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STA  TEMENT  OF  THE  PROBLEM  STUDIED 


The  goal  of  the  present  research  was  to  address  key  issues  of  vorticity  dynamics  prior 
to,  during,  and  after  the  initiation  of  dynamic  stall  and  to  develop  a  fundamental 
understanding  of  the  physical  phenomena  taking  place.  The  two  major  thrusts  of  the 
research  were  (1)  to  investigate  the  applicability  of  indicial  theory  to  the  theoretical 
prediction  of  dynamic  stall,  and  (2)  to  analyze  the  Reynolds  number  scaling  of  key  events  in 
the  dynamic  stall  process,  at  laminar  Reynolds  numbers  significantly  higher  than  could 
previously  be  computed. 

The  basic  idea  behind  the  first  portion  of  this  work  is  founded  on  recent  evidence  that 
unsteady  separated  flow  may  be  linearizable  with  respect  to  its  boundary  condition,  in  a 
manner  similar  to  linearization  techniques  used  in  unsteady  transonic  flow.  This  similarity 
(i.e.,  a  nonlinear  mean  with  time-linearizable  unsteady  component)  suggests  the 
applicability  of  powerful  analytical  tools  such  as  linear  and  nonlinear  indicial  theory.  The 
applicability  of  indicial  theory  to  key  aspects  of  the  vorticity  dynamics  of  dynamic  stall 
could  provide  the  basis  for  a  fast  prediction  method.  By  its  construction,  this  predictive 
capability  would  be  semi-analytic  and  fast-running,  making  it  well-suited  for  design  and  for 
flight  simulation  methods. 

The  second  portion  of  this  work  focused  on  the  Reynolds  number  scaling  of  incipient 
leading  edge  stall.  It  has  been  suggested  that  eruptive  plumes  of  vorticity  might  play  a 
critical  role  in  the  physics  of  vortex  initiation/ formation  during  leading  edge  dynamic  stall. 
Numerical  simulations  at  low  Reynolds  number  do  not  seem  to  adequately  predict  this 
phenomenon.  To  explore  the  possibility  of  a  "bifurcation"  in  Reynolds  number,  the  scaling 
of  incipient  laminar  separation,  of  vortex  formation,  and  of  shedding  with  respect  to 
Reynolds  number  was  investigated  in  detail .  Navier-Stokes  simulations  were  used  to  study 
a  model  problem  in  which  a  two-dimensional  airfoil  remains  stationary  at  angle-of-attack, 
but  for  which  the  leading  edge  flow  separates  as  a  result  of  an  impulsively  applied  no-slip 
boundary  condition. 


SUMMARY  OF  MOST  IMPORTANT  RESULTS 


The  research  performed  under  the  present  contract  followed  two  paths.  The  first  phase 
of  the  work  (April  92  to  August  93)  was  concerned  with  the  validation  of  indicial  theoretical 
concepts  as  applied  to  the  prediction  of  the  vorticity  dynamics  of  dynamic  stall.  The  second 
phase  of  the  work  (January  94  to  June  94)  was  an  investigation  into  high  Reynolds  number 
(up  to  800,000,  based  on  chord)  leading  edge  vortical  laminar  flow  phenomena  related  to 
stall,  based  on  a  model  problem.  The  technical  results  are  summarized  in  the  appendices. 
Appendices  A  and  B  are  copies  of  technical  papers  from  conference  proceedings. 
Appendices  C  and  D  contain  previously  unreported  results  . 
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Indicial  Theory  and  the  Prediction  of  Dynamic  Stall 

FoUowing  the  proposal  for  this  work  (Ref.  1),  the  first  portion  of  this  stt^y  began  ^th 
the  analysis  of  a  pitching  airfoil  at  a  chord  Reynolds  number  of  1,000  This  Reynolds 
number  was  found  to  be  too  low,  for  several  reasons.  The  first  one  is  that  this  case  does  not 
simulate  the  leading  edge  type  of  stall  that  is  of  practical  interest.  The  second  one  is  that  it 
does  not  exhibit  the  presence  of  eruptive  plumes  of  vorticity  which  could  be  ‘^"^jcal  to  the 
initiation  of  dynamic  stall.  The  third  reason  is  that,  at  low  Reynolds  numbers, 
flow  is  always  separated  to  varying  degrees  for  fixed  angle  of  attack,  ‘  , 

whole  idea  of  applying  indicial  theory  up  to  and  in  the  vianity  of  an  effecfave  bifurcation  of 
the  leading  edge  vbrticity  dynamics.  Here,  "bifurcation"  refers  to  a  hypothetical  s^aration 
between  a  boundary  layer-type  of  vorticity  transport  and  a  new  mode  of  vortiaty  taa^sport 
(normal  to  the  wall)  associated  with  the  appearance  of  an  eruptive  plume  of  vortiaty. 

Successful  validation  of  indicial  theory  was  immediately  achieved  at  Re  =  1,000  under 
conditions  of  high  frequency  (a)c/2U„  -  10)  and  low  amplitudes  (i.e.,  fractions  of  a  de^ee). 
However,  for  ifrger  Lplitude  motions  (Aa  =  0(5 »))  and/or  higher  Reynolds  number 
(Re  =  12,000),  appUcation  of  Unear  indicial  theory  in  its  onginal  form  failed.^The  reasons  for 
the  failure  were  as  follows;  (1)  the  flow  response  is  statically  nonlinear  when  the  airfoil  is 
beyond  static  stall;  (2)  the  incremental  asymptotic  flow  response  to  an  elemental  ^tenge  m 
the  boundary  conditions  is  inherently  unsteady  (chaotic)  at  the  higher  Reynolds  numbers. 
Indicial  theory  does  not  address  either  of  these  situations. 

In  view  of  these  findings,  the  second  portion  of  the  research  concentrated  on  the 
development  and  vaUdation  of  a  modified  version  of  indicial  theory  which  was  d^ign^^  t 
circumTent  the  above  problems.  The  resulting  characteristic/indicial  approach  uses  the 
responle^ra  reference  large  amplitude  motion  which  is  analyzed  and  processed  in  he 
Laplace  domain  to  infer  an  equivalent,  re-scaled  "indicial  f 

"characteristic  response"  of  the  system.  The  characteristic  response  can  then  be  convoluted 
with  the  time  derivative  of  the  motion  history  to  predict  the  flow  response  to  arbitrary 
schedules  of  the  boundary  condition  (a,  in  this  case).  In  essence,  successMl  Fediction  ^mg 
the  Laplace  domain  characteristic  method  signifies,  in  a  system  exhAitmg  quasi-static 
amplitude  nonlinearity,  that  the  phase  dynamics  can  still  be  treated  ^  ° 

course,  the  essence  of  the  generalized  superposition  integral 

formulation.  Alternatively,  the  new  method  can  be  interpreted  as  the  ®b-a/gbt-forward 
application  of  indicial  concepts,  using  a  re-scaled  indicial  f^^chon 
knowledge  of  the  finite  (rather  than  infinitesimal)  response  of  the  flow.  This  stretchi  g 
takes  into  account  quasi-static  nonlinearity  and  retains  the  linear  formalism  of  indic 

theory. 

By  extending  an  essentially  Unear  tool  well  into  the  nonlinear  regime,  loss  of  generality 
has  been  incurred.  In  the  present  case,  the  prediction  method  was  shown  to  be  reasonab  y 
accurate  for  a  class  of  airfoil  motions  with  arbitrary  schedule  a(t)  all  having  the  same  initial 
and  final  angles  of  attack.  Specifically,  the  time  history  of  the  vorticity  accumulation 
(defined  over  a  small  control  volume  at  the  leading  edge  of  the  airfoil) 
predicted  accurately  by  the  theory  for  different  pitch  rates  (including  oscillatory  motion), 
pr^ided  that  the  indicial  function  is  inferred  from  the  difference  between  the  flow 
responses  associated  with  any  two  airfoil  motions  a(t)  mentioned  above.  This  new 
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approach  was  validated  at  Re  =  12,000  for  a  variety  of  pitching  motions  between  a  =  10® 
and  a  =  15®. 

Under  these  conditions,  reasonable  predictive  accuracy  can  be  obtained  (on  the  order 
of  10  to  15%  on  the  vorticity  fluxes).  Typical  results  are  given  in  Appendix  A.  To  obtain 
higher  accuracy,  two  avenues  of  research  were  considered.  The  first  one  is  to  extend  the 
Laplace-domain  formalism  to  nonlinear  indicial  theory.  Efforts  in  this  direction  were 
undertaken  towards  the  end  of  the  second  funding  increment  of  this  contract.  As  a  first 
step,  these  efforts  were  specifically  directed  towards  a  piecewise  linear  version  of  the 
method.  However,  this  method  was  not  sufficiently  validated,  due  to  lack  of  time  and  to 
changes  in  research  priorities  for  the  second  phase  of  this  work  (Jan.  94  through  June  94). 
The  second  avenue  of  research  involved  the  separation  of  deterministic  (motion-related) 
and  chaotic  portions  of  the  predicted  signals.  One  of  the  difficulties  associated  with  the 
application  of  the  prediction  method  to  laminar  single-realization  cases  was  the  presence  of 
significant  chaotic  fluctuations  associated  with  vorticity  dynamics  downstream  of  the 
leading  edge  region.  The  chaotic  unsteady  effects  are,  to  a  large  degree,  uncorrelated  with 
the  change  in  boundary  condition,  Aa(t).  In  particular,  the  method  can  only  predict  the 
deterministic  component  of  the  flow  response.  For  instance,  for  sufficiently  large  times,  the 
predicted  increment  response  A/i  (between  two  finite-amplitude  motions  having  same  final 
angle-of-attack)  equals  zero,  since  the  input  Aa  to  the  method  is  zero  for  a  time  exceeding 
the  length  of  the  convolution.  Yet,  for  sufficiently  large  Reynolds  number,  A/i(t)  ^  0.  This  is 
expected  of  deterministic  methods  with  finite  "memory."  A  similar  problem  was  found  in 
the  artificial  neural  network  predictions  of  Faller  et  al.  (Ref.  2),  for  large  times.  Since  this 
aspect  of  the  prediction  cannot  be  circumvented,  a  search  for  candidate  dynamic  stall 
"indicators"  was  undertaken.  The  basic  idea  is  that,  if  A//(t)  is  small,  the  prediction  error  will 
also  be  small.  Thus,  a  suitable  indicator  must  have  the  property  that,  in  response  to  the 
airfoil  motion,  its  time-history  should  exhibit  a  strong  degree  of  determinism  with  respect 
to  Aa(t). 

The  applicability  of  indicial  theory  hinges  on  the  existence  of  flow  quantities  whose 
temporal  behavior  remains  smooth  throughout  dynamic  stall.  A  result  of  this  research  is 
that  vorticity  accumulation  appears  to  be  a  suitable  candidate  as  a  dynamic  stall  "indicator" 
in  which  the  deterministic  (motion-related)  component  is  dominant.  Vorticity  accumulation 
also  happens  to  be  predicted  well  by  indicial  theory  as  modified  in  this  work.  The 
applicability  of  the  modified  theory  indicates  that  the  phase  dynamics  are  linear,  and  that  a 
generalized  superposition  integral  can  be  used  to  construct  phase  lags  from  elementary 
responses. 

Reynolds  Number  Effects  on  a  Model  Problem  of  Leading  Edge  Stall 

The  study  of  a  model  problem  was  considered  in  which  an  airfoil  is  initially 
maintained  at  a  fixed  angle-of-attack  (a  =  15®).  A  converged  Navier-Stokes  solution  which 
is  free  of  separation  is  obtained  by  using  the  zero-vorticity  slip  boundary  condition.  At 
some  time  t  =  0,  the  no-slip  boundary  condition  is  impulsively  turned  on.  The  resulting 
developing  flowfield  is  characterized  by  the  complete  absence  of  trailing  edge  stall,  the 
existence  of  a  well-defined  time  origin,  and  the  independence  of  the  initial  condition  on 
Reynolds  number.  This  problem  was  originally  conceived  of  in  the  initial  part  of  the  work, 
in  order  to  mimic  certain  physical  aspects  of  high  Reynolds  number  flow  in  low  Reynolds 
number  simulations.  In  particular,  this  approach  circumvents  the  "contamination"  of  the 
leading  edge  flow  due  to  reverse  flow  emanating  from  the  trailing  edge  region. 
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Incidentally,  the  properties  of  the  flowfield  associated  with  this  model  problem  are  also 
advantageous  for  the  detailed  study  of  the  very  initial  stages  of  dynamic  stall,  as  well  as  for 
the  investigation  of  scaling  effects  with  Reynolds  number. 

High  Reynolds  number  laminar  flow  remains  a  difficult  area  of  research.  From  a 
computational  point  of  view,  the  currently  accepted  highest  chord  Reynolds  number  at 
which  grid-independent  results  can  be  obtained  is  Re  =  100,000.^  With  the  present  model 
problem  and  using  a  high-order  differencing  code,  grid-converged  solutions  of  the  two- 
dimensional  Navier-Stokes  equations  were  obtained  up  to  Re  =  800,000.  The  importance  of 
the  high  Reynolds  number  laminar  research  for  dynamic  stall  stems  from  high  Reynolds 
number  asymptotics^  and  interacting  boundary  layer  research^'^  which  suggest  that 
eruptive  plumes  of  vorticity  may  form  in  the  leading  edge  boundary  layer  and  act  as  a 
triggering  mechanism  leading  to  the  formation  of  the  primary  dynamic  stall  vortex.  Such 
eruptions  (which  differ  from  the  standard  vortex-induced  eruption)  have  not  yet  be^  foimd 
in  numerical  solutions  of  the  Navier-Stokes  equations.  A  prevalent  point  of  view^  (and  a 
strong  possibility)  is  that,  in  an  Eulerian  reference  frame,  the  mesh  spacing  to  resolve  this 
phenomenon  cannot  be  afforded.  Another  possibility  is  that  some  form  of  bifurcation  takes 
place  at  some  intermediate,  laminar  Reynolds  number.  Since  most  numerical  calculations 
are  at  low  Reynolds  number,  the  effects  of  viscosity  may  dissipate  any  tendency  to 
spontaneously  form  these  eruptive  plumes.  The  main  idea  behind  an  analysis  of  the 
Reynolds  scaling  of  key  events  in  the  dynamic  stall  process  is  that  the  emergence  of  new 
physics  -even  if  poorly  resolved-  would  result  in  modified  scaling  properties  with  respect  to 
Reynolds  number. 

A  high-order  accuracy  finite  difference  code  was  used  for  this  part  of  the  investigation. 
The  code  uses  eighth-order  spatial  accuracy  on  the  convective  fluxes,  on  the  metrics,  and  on 
the  artificial  dissipation,  while  retaining  fourth-order  accuracy  on  the  viscous  terms.  The 
numerical  method  allows  for  the  accurate  resolution  of  a  wavenumber  range  in  each  spatial 
direction  equivalent  to  80%  of  the  wavenumber  range  for  pseudo-spectral  methods.  Thus 
the  required  number  of  grid  points  is  substantially  reduced,  compared  to  a  standard  second- 
order  code.  A  more  detailed  description  of  the  code  can  be  found  in  Appendices  C  and  D. 
High  order  accuracy  was  a  critical  element  which  permitted  the  extension  of  the  model 
problem  studies  to  Reynolds  numbers  considerably  larger  than  previously  possible. 

At  least  two  distinct  stages  of  the  incipient  separation  process  were  found.  Stage  I 
corresponds  to  the  formation  of  a  separation  zone  ("lifting"  of  the  thin  primary  vorticity 
layer)  which  is  initially  symmetric  about  the  location  of  maximum  botmdary  layer  thickness. 
This  stage  was  found  to  exhibit  qualitatively  self-similar  behavior  over  the  investigated 
Reynolds  number  range  (50,000  <  Re  <  800,000,  M„  =  0.2),  and  ends  with  the  accumulation 
of  vorticity  in  the  reattachment  region  of  the  separated  vorticity  layer.  In  Stage  II,  the 
accumulated  vorticity  develops  into  a  vortex  which  subsequently  induces  a  sequence  of 
rapid  events,  including  the  formation  of  secondary  flow  and,  at  high  Reynolds  number, 
successive  shear  layer  vortices  which  themselves  induce  secondary  and  tertiary  flow 
patterns.  This  stage  is  the  most  difficult  to  compute,  and  issues  of  adequacy  of  the 
numerical  resolution,  diffusion  of  the  vortices,  and  steepness  of  the  induced  pressure 
gradients  are  critical.  At  high  Reynolds  numbers,  a  second  order  method  would  be 
extremely  expensive  if  it  were  to  capture  the  essential  physics. 

At  the  beginning  of  this  portion  of  the  research,  a  number  of  important  research  issues 
were  raised,  in  particular  with  respect  to  the  relevance  and  applicability  of  high  Reynolds 
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number  asymptotic  theory  and  interacting  boundary  layer  calculation  results  (see,  e.g.. 
Ref.  5).  The  theory  suggests  the  development  of  a  tightly  focused  narrowing  band  eruption 
in  the  boundary  layer.  This  eruption  is  a  generic  boundary  layer  response  to  a  moving 
adverse  pressure  gradient  zone  outside  the  boundary  layer.  The  presence  of  a  vortex- 
induced  eruption  in  Stage  11  is  well-established,  and  has  been  observed  in  the  present 
calculations  at  Re  =  100,000, 200,000, 400,000,  and  800,000.  Knight  and  Choudhuri's  results^ 
exhibit  similar  behavior  at  a  Reynolds  number  as  low  as  Re  =  10,000.  There  is  also  reason  to 
believe  that  this  type  of  eruption  would  presumably  also  have  taken  place  here  at 
Re  =  50,000,  had  the  calculations  been  carried  out  further  in  time.  However,  no  eruption 
was  ever  observed  during  Stage  I  of  the  flow  development,  nor  could  an  eruption 
(bifurcation)  be  inferred  from  the  basic  Reynolds  number  scaling  properties  of  the  flow. 

The  present  results  indicate  a  smooth  scaling  behavior  over  the  Reynolds  number 
range  50,000  s  Re  <  800,000  in  Stage  1,  and  the  well-defined  emergence,  from  the  lifting 
vorticity  layer,  of  a  topological  "center"  which  precedes  vortex  formation.  Qualitatively 
similar  behavior  was  observed  at  all  investigated  Reynolds  numbers.  Furthermore,  when 
taking  into  account  the  viscous  scaling  of  the  Stokes  layer  and/ or  boundary  layer  upstream, 
the  specific  dynamical  events  associated  with  incipient  stall  (such  as  the  times  of  onset  of 
primary  and  secondary  separation  bubbles,  and  of  vortex  formation)  were  found  to  scale 
according  to  simple  power  laws  of  the  Reynolds  number.  Re"'*,  with  0.11  </}  ^  0.45.  Self¬ 
similarity  was  also  observed  for  the  time  of  occurrence  of  vortex-induced  eruptions 
(marking  the  onset  of  shedding  of  the  primary  stall  vortex),  between  the  Reynolds  numbers 
of  50,000  and  400,000.  At  Re  =  800,000,  new  physics  associated  with  complex  interactions 
between  vortices  take  place  in  the  vicinity  of  the  wall  before  the  time  of  eruption.  However, 
no  qualitative  changes  were  observed  in  Stage  I. 

The  present  results  contradict  the  hypothesis  that  a  form  of  bifurcation  must  take  place 
at  some  intermediate  laminar  Reynolds  number  during  Stage  1,  i.e.,  before  the  formation  of 
a  vortex.  To  reconcile  these  results  with  the  implications  of  interacting  boundary  layer 
theory,^  one  of  several  conjectures  may  be  made.  These  conjectures  are  listed  below. 

The  first  one  concerns  the  Reynolds  number  of  the  simulations.  Although  Re  =  800,000 
goes  beyond  the  traditionally  accepted  Re  =  100,000  limit  for  laminar  calculations  of 
dynamic  stall,  the  fact  remains  that  this  is  a  low  Reynolds  number  flow.  The  Reynolds 
number  range  based  on  leading  edge  radius,  b,  was  3,250  <  Re^,  <  12,800.  Based  on 
boundary  layer  thickness  (S^l^)  and  freestream  velocity,  the  Reynolds  number  range  was: 
46  <  Rej  <  186.  These  are,  of  course,  extremely  low  Reynolds  numbers,  compared  to  those 
considered  in  the  high  Reynolds  number  asymptotics  and  interacting  boundary  layer 
studies. 

The  second  conjecture,  also  related  to  the  Reynolds  number,  is  that  if  a  Reynolds 
number  bifurcation  in  Stage  1  should  occur,  it  will  take  place  at  a  Reynolds  number  higher 
than  those  investigated  in  this  study.  It  is,  therefore,  likely  that  transition  and  turbulence 
will  be  present.  The  effects  of  transition  and  turbulence  could  overwhelm  the  bifurcation 
phenomenon,  thus  making  it  very  difficult  to  discern.  Calculations  at  chord  Reynolds 
numbers  larger  than  800,000  were  not  attempted  in  this  study,  in  part  because  transition  to 
turbulence  may  take  place  within  the  thin  free  shear  layers.  Although  the  present  resolutioii 
is  sufficient  to  resolve  the  primary  Kelvin-Helmholz  instability  (at  least,  for  the  primary  free 
shear  layer),  the  present  calculations  are  no  longer  accurate  when  three-dimensional 
instabilities  are  present. 
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The  third  conjecture  is  that  a  self-focusing  plume  of  vorticity  must  exist  at  the  leading 
edge  within  the  investigated  Reynolds  number  range.  In  this  scenario,  the  reason  why  this 
was  not  observed  is  because  the  phenomenon  is  subscale  and  cannot  be  resolved  by  Eulerian 
means.  Although  this  is  a  theoretical  possibility,  experience  teaches  us  that  some  sign  of 
"activity,"  be  it  instability  or  other  phenomenon,  would  have  been  visible  in  the  calculations, 
especially  since  an  eighth-order  code  was  used.  No  such  sign  was  discovered  in  any  of  the 
present  calculations.  Grid  convergence  was  attained  at  all  Reynolds  numbers,  and  further 
grid  refinements  were  not  attempted.  The  prospect  of  further  increasing  grid  resolution  to 
repeat  the  same  results  was  considered  unproductive  and  risky;  therefore,  this  avenue  of 
research  was  not  pursued. 

The  fourth  conjecture  pertains  to  the  absence  of  a  moving  region  of  adverse  pressure 
gradient.  The  physics  of  formation  of  eruptive  plumes  of  vorticity  are  described  by  the  high 
Reynolds  number  asymptotics  and  interacting  boundary  layer  studies.'*'  These  studies 
specifically  address  the  response  of  the  boundary  layer  to  a  moving  zone  of  adverse 
pressure  gradient.  In  the  present  semi-impulsive  problem  at  fixed  angle  of  attack,  there  is 
very  little  temporal  change  in  the  leading  edge  pressure  gradient.  Except  for  an  initial 
relaxation  of  the  pressure  gradient  after  imposing  the  no-slip  condition,  the  region  of 
adverse  pressure  gradient  does  not  undergo  any  significant  motion.  This  is  different  from  a 
pitching  airfoil  situation.  Whether  or  not  this  is  significant  is  something  that  could  be 
systematically  investigated  by  superimposing  a  pitching  motion  to  the  present  problem. 
However,  this  route  was  not  pursued  within  the  time  frame  of  this  study,  for  two  reasons. 
The  first  one  is  the  fact  that  vortex-induced  eruptions  take  place,  even  though  the  motion  of 
the  vortex  is  very  slight:  the  temporal  strengthening  of  the  circulation  associated  with  the 
vortex  appears  to  be  by  far  the  largest  contributor.  The  second  reason  is  that,  if  one 
considers  the  time  scales  involved  for  the  entire  process  to  complete  (say,  up  to  the  point  of 
shedding  of  the  primary  stall  vortex),  these  correspond  to  only  a  fraction  of  a  degree  for  a 
realistic  airfoil  motion.  For  example,  for  a  hypothetical  oscillatory  motion  at  Re  =  800,000  at 
a  reduced  frequency  k  =  wc/2U  =  0.02  with  a  motion  amplitude  of  5 ",  the  change  in  angle  of 
attack,  Ao,  from  initiation  to  shedding,  equals  0.014®.  This  would  not  be  expected  to 
produce  any  significant  motion  of  the  zone  of  adverse  pressure  gradient. 

Finally,  the  results  presented  here  suggest  that  the  behavior  of  the  flow  started  to 
change  at  Re  =  800,000.  Although  these  differences  occurred  in  Stage  11,  it  is  useful  to 
consider  the  broader  picture  and  attempt  to  speculate  as  to  what  happens  at  higher 
Reynolds  numbers.  At  Re  =  800,000,  a  vortex-induced  eruption  takes  place  virtually 
underneath  a  merged  primary  vortex,  as  opposed  to  taking  place  between  consecutive 
primary  and  shear  layer  vortices.  The  result  is  a  fundamental  change  from  a  bulk 
convection  of  the  vortex  system  at  low  Reynolds  numbers  to  a  sudden  break-away  from  the 
surface  at  a  point  in  time  when  the  entire  vortex  system  is  confined  to  the  leading  edge 
region.  As  the  Reynolds  number  becomes  larger,  the  artificial  distinction  between  Stages  1 
and  II  will  be  blurred  and  not  realistic.  Transition  effects  will  likely  replace  the  elongated 
laminar  bubble  of  Stage  I  by  transitional  and  turbulent  flow  structures  in  the  close  proximity 
of  the  wall.  If  this  is  the  case,  the  broader  issue  is  to  determine  the  conditions  under  which 
eruption  of  vorticity  occurs  if,  as  suggested  by  the  present  results,  such  an  eruption  acts  as 
the  trigger  for  vortex  shedding. 
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TOWARDS  A  SEMI-ANALVnC  TOOL  AIAA  94-0537 

FOR  THE  PREDICTION  OF  DYNAMIC  STALL 


Patrick  H.  Reisenthel 
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ABSTRACT 

Direct  numerical  simidations  are  used  to  analyze 
in  detail  the  vorticity  d)mamics  of  the  leading  edge 
region  of  a  NACA0012  airfoil  pitched  about  its 
1  /4  chord.  The  results  presented  in  this  paper  illustrate 
how  the  formalism  of  indicial  theory  can  be  used  to 
predict  the  integrated  vorticity  fluxes  and  the  vorticity 
acciunulation  during  unsteady  maneuver.  In  particular, 
the  flow  response  to  large  amplitude  non-linear  motions 
is  shown  to  be  predicted  reasonably  accurately,  provided 
that  the  indicial  functions  of  the  flow  are  inferred  in  the 
Laplace  domain  and  stretched  to  account  for  quasi-static 
non-linearity.  The  implication  of  this  work  is  the 
possibility  of  developing  a  fest  semi-analytical  prediction 
method  for  incipient  leading  edge  stall,  which  will  be 
accurate  within  certain  classes  of  maneuvers. 


NOMENCLATURE 

Symbols 

A  Amplitude  transfer  function  (A  =  A(«)) 

Afl  Quasi-static  amplitude  (Aq  =  A(0)) 

c  Airfoil  chord 

amst.  Constant 

CFD  Computational  Fluid  Dynamics 
DNS  Direct  Numerical  Simulation 
/  Flow  response 
F  Frequency 

k  Reduced  frequency  (J:  =  (2jrFc)/(2UJ) 

L  Laplace  transform 

L. E.  Leading  Edge 

M, .  Freestream  Mach  number 

Re  Chord  Reynolds  number  (Re  =  U„c/vJ 
s  Curvilinear  (tangential)  airfoil  coordinate 
(0.0  sss  1.0) 

s  Complex  Laplace  variable  ^ 

t  Non-dimensional  time  (t  =  t  c) 

T  Ramp  time  (non-dimensionalized  with  U,. 
andc) 

u  Streamwise  velocity 
U„  Freestream  velocity 
V  Normal  velocity 
X  Streamwise  coordinate 

X  Generic  spatial  coordinate 
y  Normal  coordinate 

yj^  Coordinate  locally  perpendicular  to  the  airfoil 
surface 

a  Angle  of  attack  (in  degrees  unless  otherwise 
specified) 


£  Forcing  function 

A  Variation 

A#  Phase  transfer  function  (A^  =  A^fo))) 

•  Total  vorticity  flux  (convective  and  viscous) 

a  Vorticity  flux  angle 

y  Flow  response  (generic  notation) 

ft  g  Indicial  function  (generic) 

V  Fluid  kinematic  viscosity 

p  Fluid  density 

T  Time  lag 

a  Reduced  atjgular  frequency  («=2*Fc/UJ 

n  Vorticity 

f  Time  variable  for  a  predeflned  motion  o(f) 

Stibsaipts 

e  Evaluated  at  the  outside  edge  of  the  botmdary 

layer 

o  Initial  condition 

n  Normal 

t  .  Tangential 
<>>  Freestream  value 

Superscripts 

-1  Inverse 

-  Indicial  function 

*  Accumulation 

+  Non-dimensional 

''  Dimensional 


INTRODUCTION 

The  past  two  decades  have  seen  a  considerable 
amount  of  resources  spent  towards  resolving  the 
problem  of  dynamic  stall.  Despite  some  progress  in  both 
numerical  prediction^  and  experimental^^  techniques, 
there  are  still  too  many  unknowns  to  ensiue  successful 
transfer  of  this  knowledge  from  basic  research  to  appli^ 
research  and  exploratory  development.  A  case  in  point 
is  the  design  of  helicopter  rotor  blades,  which  makes 
little  use  of  dynamic  stall  results.  Ideally,  a  successful 
blade  design  might  be  one  which  could  benefit  from  the 
dynamic  lift  effect,  without  experiencing  the  enlunced 
torsional  fatigue  caused  by  shedding  of  a  d5mamic  stall 
vortex. 

The  iiuibility  of  the  research  community  to  provide 
designers  involved  in  the  aerod)mamics  of  unsteady 
lifting  surfaces  with  the  appropriate  tools  is,  without  a 
doubt,  related  to  the  difficulty  of  the  problem  at  hand 
and  has  been  a  source  of  constant  frustration.  The 
current  understanding  of  the  complex  physics  involved 
in  dynamic  stall  is  limited  by  such  difficult  and  inter- 
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related  issues  as  three-dimensionality^  and  transition  to 
turbulence  (e.g.,  moving  wall®  effects),  as  well  as 
fimdamental  scaling  issues  with  Reynolds  number®  and 
additional  physics  introduced  by  compressibility  effects.^ 

The  basic  idea  behind  the  work  presented  in  this 
paper  is  founded  on  recent  evidence^  that  unsteady 
separated  flow  may  be  linearizable  with  respect  to  its 
boundary  condition,  in  a  manner  similar  to  unsteady 
transonic  flow.®  This  similarity  (i.e.,  a  non-linear  mean 
flow  with  time-linearizable  unsteady  component) 
suggests  the  applicability  of  powerful  analytical  tools 
such  as  linear  and  non-linear  indicial  theory.*'^®''*  The 
applicability  of  indicial  theory  to  d3mamic  stall,  or  its 
applicability  to  some  key  aspect  of  the  vorticity 
dynamics  associated  with  the  dynamic  stall  process, 
could  provide  the  basic  element  of  a  fast  prediction 
method.  By  its  construction,  this  predictive  capability 
would  be  semi-analytic,^'®  making  it  well-suited  for 
design  and  for  flight  simulation  techniques. 

OBJECTIVE  AND  APPROACH 

The  objective  of  this  work  is  to  attempt  to  either 
prove  or  disprove  the  applicability  of  the  indicial 
theoretical  concept,  as  applied  to  the  dynamic  stall 
process.  The  approach  is  to  use  direct  numerical 
simulations  (DNS)  of  dynamic  stall  for  a  well- 
documented  two-dimensional  configuration:  that  of  a 
NACA0012  airfoil  in  pitch  about  the  1/4  chord  position. 
Since  indicial  theory  (summarized  below)  is  applied  here 
to  the  physics  of  incipient  unsteady  separation,  it  is 
essential  that  turbulence  models  (which  would 
unavoidably  "contaminate"  the  physics  which  are  the 
very  object  of  this  analysis)  be  eliminated.  Therefore,  as 
a  first  step,  the  Reynolds  numbers  considered  here  are 
low  (Re  =  1,000,  Re  =  12,000,  and,  in  a  model  problem. 
Re  =  100,000),  corresponding  to  laminar  flow  corulitions. 

The  analysis  is  based  on  the  hypothesis  that  the 
cause  of  dynamic  stall  is  linked  to  the  disparate 
characteristic  times  for  the  various  physical  processes 
that  take  place.  The  effect  of  the  disparity  in  time  scales 
is  such  that,  for  a  given  pitch  rate,  the  following 
processes  occur:  i)  at  some  critical  time,  the  vorticity 
accumulated  throughout  the  airfoil  motion  exceeds  the 
maximum  allowable  for  vorticity  transport  through  an 
attached  boundary  layer  model;  at  that  time,  normal 
vorticity  transport  component  becomes  important,  and 
ii)  the  normal  vorticity  transport  compwnent  is  confined 
to  a  vortex  feeding  sheet  or  a  sharply  defined  shear  layer 
which  may  lead  to  the  formation  of  the  vortex.  In  steady 
flow  the  separation  may  start  from  a  shear  layer  but 
rapidly  becomes  chaotic.  In  the  dynamic  case  the  time 
lags  may  allow  the  sustained  development,  or  stability, 
of  the  vortex  sheet,  such  that  a  vortex  is  formed.  Once 
formed,  the  vortex  may  interact  with  the  boundary  layer 
below  it,  eventually  giving  rise  to  an  eruptive  plume  of 
vorticity.  The  result  of  this  eruption  is  to  free  the 
dynamic  stall  vortex  from  its  feeding  sheet,  which  is 
followed  by  the  rapid  convection  of  the  vortex  over  the 
surface. 


A  key  addition  to  this  idea  comes  from  recent 
work^  in  which  it  was  found  that  if  an  impulsive  change 
is  given  to  the  flow  around  a  cylinder,  the  transient 
value  of  the  shed  vorticity  flux  at  short  times  is 
considerably  larger  than  the  value  for  steady  flow.  It  is 
noteworthy  that  the  hypothesis  that  time  lags  are 
responsible  for  the  delayed  breakup  of  the  vortex  sheet  is 
consistent  with  the  findings  of  Henk  et  al.*^  who 
experimentally  examined  imsteady  separated  flows  in  a 
water  tunnel. 


Unlike  pure  Computational  Fluid  Dynamical 
(CFD)  approaches,  which,  like  experiments,  produce 
numbers  (but  not  necessarily  insight),  and  are  too 
expensive,  the  semi-analytical  approach  outlined  below 
provides  a  guide  to  physical  understanding  since  it 
allows  the  identification  of  cause  and  effect  relationships. 


Indicial  method.  The  present  approach  is  ba^ed  on  the 
concept®  that  a  characteristic  flow  variable  f(x,t),  which 
describes  the  state  of  the  flow,  can  be  linearized  with 
r^pect  to  its  boundary  condition,  c(t),  if  the  variation  of 
f(x,t)  is  a  smooth  function  of  £(t).  Since,  in  general, 
separation  is  a  non-linear  process,  the  function  f(x,t) 
must  be  such  that  its  variation  with  c(t)  is  smooth 
through  separation;  this  allows  the  representation  of 
f(x,t)  in  a  Taylor  series  about  some  value  of  c,  thus 

+  .  .  .  (1) 

If  f(x,0)  is  zero  (a  zero  initial  condition),  then  an 
approximate  solution  is 

(2) 

e=eo 

Equation  (2)  is  an  approximation  which  becomes  more 
accurate  as  A£  0.  Also,  Eq.  (2)  is  exact  if  f(x,t)  is  a  linear 
function  of  c(t).  If  (2)  is  sajsfied,  then  it  can  be  shown® 
that  a  formal  solution  for  f(x,t)  is 

t 

f(i.t)  =  ?(i,t)e{0)  -  j  ?(x.T)dT  (3) 

where 


£{i,t)  =  A*  ^ 


f(x.t)  =  f(x,0)  +  Ac  ^ 


?(x,t) 


it 

le 


e=So 


(4) 


Hence,  if  the  forcing  function  (i.e.,  the  boundary 
condition  e)  is  known  and  if?  (the  "indicial  function")  is 
known  from  some  computation,  then  Eq.  (3)  gives  the 
value  of  f(x,t)  for  any  schedule  of  boundary  conditions  c(t) 
without  the  need  to  compute /from  first  principles.  This 
has  the  effect  of  reducing  computational  reqviirements 
considerably.  FurAermore,  Eq.  (3)  is  a  semi-analytic 
relation  between  f(x,t)  and  its  boundary  condition  e(t) 
which  can  give  insight  into  the  physical  process 
involved. 


If  applicable,  the  formalism  of  the  theory 
represented  by  Eq.  (3)  allows  insight  into  the  character  of 
dynamic  stall  even  though  the  actual  details  of  the 
response  are  unknown.  The  analysis  outlined  above 
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illustrates  the  power  of  the  method,  a  power  that  is 
greatly  enhanced  if  details  of  the  response  are  computed. 

It  may  be  difficult  to  comprehend  how  a  linear 
equation  can  represent  the  non-linear  behavior  of 
separation,  and  therefore  it  is  useful  to  illustrate  how  this 
concept  represents  an  example  of  unsteady  separated 
flow;  in  this  case,  the  flow  around  an  oscillating  cylinder. 
This  example  is  based  on  work  published  in  Ref.  7  which 
shows  that  it  is  possible  to  predict  the  time-dependent 
characteristics  of  separation  location  and  vorticity  flux, 
based  on  the  knowledge  of  the  impulse  response  of  the 
flowfield  only.  Previous  experience  with  the  appliration 
of  indicial  theory  to  the  prediction  of  non-linear 
transonic  flow^  indicates  that  whUe  the  overall  process 
may  be  non-linear,  quantities  which  vary  smoothly 
throughout  the  change  in  boundary  condition  can  be 
accurately  predicted  using  indicial  theory.  This  applies 
in  particular  to  integral  quantities  such  as  overall  loads, 
provided  that  the  topology  of  the  flowfield  does  not 
change  (e.g.,  appearance  or  disappearance  of  shock 
waves).'^ 

Hence,  the  applicability  of  the  indicial  method  to 
dynamic  stall  hinges  on  the  existence  of  quantities  which 
remain  smooth  throughout  the  phenomena  associated 
with  dynamic  stall.  As  will  be  shown  in  the  Results 
Section  below,  vorticity  fluxes  in  the  leading  edge  region 
appear  to  satisfy  this  requirement.  Furthermore,  even  in 
the  hypothetical  case  (not  treated  here)  where  a  flow 
quantity  of  interest  should  undergo  some  form  of 
bifurcation,  the  only  requirement  for  predicting  the  onset 
of  dynamic  stall  is  that  its  unsteady  component  behave 
piecewise  linearly  on  one  side  of  (i.e.,  prior  to)  this 
bifurcation.  For  example,  if  the  so-called  vorticity 
eruption  phenomenon^^'^®  at  the  leading  edge 
constitutes  a  bifurcation  point,  then  the  previous 
theoretical  developments  apply  only  to  the  physics  of  the 
flow  up  to  the  vorticity  eruption  process. 

In  order  to  effectively  predict  dynamic  stall  using 
indicial  theory,  two  separate  but  related  objectives  must 
be  achieved.  First,  the  applicability  of  indicial  theory 
must  be  demonstrated.  Second,  and  in  the  event  that  the 
demonstration  is  successful,  proof  must  be  given  that  the 
flow  quantity(ies)  to  which  the  theory  applies  are  indeed 
precursors  of  the  phenomenon  itself.  While  some  results 
towards  this  second  objective  (not  shown  here)  are 
currently  being  obtained  as  part  of  the  contract 
supporting  this  work,  the  limited  obje^ve  of  the  pre^t 
paper  is  to  prove  or  disprove  the  applicability  of  indicial 
theory  to  the  problem  of  incipient  dynamic  stall. 


RESULTS 

All  calculations  presented  in  this  paper  have  beeii 
performed  with  a  modified  version  of  NASA  Ames' 
ARC2D^®  time-implicit  finite  difference  computer  code. 
Two  computational  grids  were  used:  a  169  x  51  grid  for 
the  Re  =  1,000  case  and  a  169  x  80  for  the  Re  =  12,000 
case.  Details  of  the  latter  grid  (O-mesh)  are  shown  in 
Fig.  1  for  the  trailing  and  leading  edge  regions  of  the 


airfoU.  Fig.  1(b)  illustrates  the  tangential  clustering  of 
the  mesh  in  the  region  between  0  and  2  percent  of  chord 
fiom  the  leading  edge.  Also  shown  (grey  shaded  area)  is 
a  control  volume  over  which  integrated  tangential  and 
normal  fluxes  of  the  vorticity  are  computed  and  stored  at 
each  time  step.  A  typical  result  (in  which  only  a  subset 
of  flow  vectors  has  been  shown  for  clarity)  appears  in 
Fig.  2. 

Although  most  of  the  quantitative  results 
presented  in  this  paper  are  for  the  Re  =  1,000  and 
Re  =  12,000  cases  (M.  =  0.2),  sample  results  for  a  model 
problem  related  to  dynamic  stall  (Re  =  100,000,  =  0.1) 

will  be  presented  first.  The  purpose  of  this  model 
problem  (described  below)  is  to  remove  some  of  the 
complications  and  simplify  the  problem,  in  order  to 
isolate  relevant  physics  and  better  understand  the 
vorticity  dynamics  associated  with  "incipient"  dynamic 
stall. 

Model  Problem 

In  the  absence  of  compressibility  effects,  vorticity 
may  only  enter  the  flowfield  by  viscous  diffusion  tom 
the  boundary  of  the  airfoil.  The  amount  of  vorticity  is 
dictated  by  the  tangential  pressure  gradient  in  the 
vicinity  of  the  leading  edge.l7,i8  jYds  pressure  gradient 
would  be  maximized  at  any  angle  of  attack  if  the  flow 
remained  completely  attached  such  as  in  potential  flow. 
In  reality,  the  pressure  gradient  is  "relieved"  by  flow 
separation,  which  acts  to  augment  the  effective  leading 
edge  radius  "seen"  by  the  outer  potential  flow.  TWs  idea 
is  also  based  on  corroborating  evidence  in  stall 
experiments  at  high  pitch  rate,  which  show  that  the 
formation  of  a  leading  edge  recirculation  is  delayed,  thus 
creating  a  flowfield  at  the  leading  edge  which  initially 
"resembles"  potential  flow. 

To  approximate  this  condition,  vorticity-free 
boundary  conditions  were  implemented,  resulting  in 
steady  initial  conditions  which  were  fi«e  of  separation 
(see  for  example  Fig.  3).  At  some  initial  time  derioted 
t  =  0,  the  no-slip  boundary  condition  is  impulsively 
applied.  The  vorticity  field  that  rapidly  develops  after 
the  imposition  of  the  no-slip  condition  results  in 
dynamic  stall  within  a  fraction  of  the  conv^five  time 
scale.  The  existence  of  a  well-defined  time  origin  for  the 
various  dynamical  processes,  the  absence  of  a  laminar 
bubble  prior  to  stall,  and  the  absence  of  any  trailing  edge 
separation  are  the  k^  ingredients  of  the  model  problem. 
In  addition,  the  initial  conditions  are  independent  of 
Reynolds  number.  The  resulting  unsteady  flowfield 
(which  will  be  the  topic  of  a  forthcoming  AIAA  Paper*’') 
shares  many  of  the  characteristics  of  high  Rejmolds 
number,  high  pitch  rate  d3mamic  stall,  with  the  exception 
of  moving-wall  and  apparent  mass  eff^. 

The  results  shown  in  Figs.  4  and  5  illustrate  the 
evolution  of  the  vorticity  field  at  successive  instants  of 
time  shortly  after  turning  on  the  no-slip  boundary 
condition.  Detailed  analyses  of  the  unsteady  flowfield  at 
a  =  45*  and  o  =  30«  (not  shown  here)  indicate  that 
leading  edge  separation  occurs  extremely  rapidly,  and 
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that  the  process  is  essentially  "smooth"  in  terms  of 
vorticity  dynamics.  The  surface  flow  topology 
corresponding  to  o  =  30*  is  shown  in  Fig.  6(a) 
(Re  =  100,000),  where  the  lines  demarcate  regions  of 
forward  and  reverse  flow  along  the  airfoil  surface.  This 
is  contrasted  to  the  Re  =  1,{X)0  case  (Fig.  6(b),  same  scale) 
which  exhibits  neither  shedding  (indicated  by  oblique 
pairs  of  lines  in  Fig.  6(a))  nor  secondary  separation  over 
the  shown  time  interval  (0.0  s  tU^c  s  0.1).  Secondary 
separation  at  Re  =  1,000  occurs  some  later  time  (see 
Fig.  7(b))  and,  hence,  the  overall  time  scale  for  the  onset 
of  dynamic  stall  is  larger  both  at  the  lower  Reynolds 
numbers  (see  Fig.  6)  and  for  lower  angles  of  attack  (note 
the  doubling  of  the  horizontal  scale  between  the  cases  in 
Fig.  7).  Such  conditions  are  better  suited  for  the  purpose 
of  detailed  temporal  analyses  of  vorticity  fluxes  and 
accumulation.  Consequently  the  remainder  of  this  paper 
is  devoted  to  the  low  Reynolds  number  flowfields 
(Re  =  1,000  and  Re  =  12,000)  as  a  test  bed  for  the 
applicability  of  indicial  theory  to  the  prediction  of 
incipient  dynamic  stall. 

A  comparison  of  the  time  evolution  of  tangential 
and  normal  fluxes  of  the  vorticity  (for  the  control  volume 
depicted  in  Fig.  1(b))  is  shown  in  Fig.  8.  Each  line 
corresponds  to  the  flux  at  a  given  (tangential  or  normal) 
position;  the  arrow  indicates  the  direction  of  increasing 
coordinate.  The  fluxes  are  seen  to  peak  and  decay 
smoothly.  An  equivalent  "flux  angle,"  0,  at  a  fixed 
tangential  location  central  to  the  control  volume  can  be 
defined  at  each  of  the  normal  coordinates  as 
0  =  tan‘*(0n/0t).  Vorticity  accumulation  (not  shown 
here)  peaks  at  approximately  tU./c  =  0.003.  This  is 
followed  by  a  gradual  "turn"  (see  Fig.  8(c))  of  the 
convective  vorticity  flux  from  a  direction  tangential  to 
the  surface  (0  =  0*)  to  an  inclination  0  approaching  45*. 
While  qualitatively  similar  peak-and-decay  behavior  of 
both  normal  and  tangential  vorticity  fluxes  was  also 
observed  for  the  case  of  real  dynamic  stall  during  steady 
pitch-up  (da*/dt  =  0.6),  the  overall  time  scales 
(consistent  with  the  results  of  Fig.  8)  are  typically  larger 
than  in  the  model  problem  at  a  =  45  • . 

Validation  of  Indicial  Theory 

Having  tentatively  established  the  relative 
smoothness  of  the  leading  edge  vorticity  dynamics  prior 
to,  during,  and  immediately  after  the  onset  of  dynamic 
stall,  we  now  turn  to  some  direct  attempts  at  the 
validation  of  indicial  theory.  This  section  is  divided  into 
three  parts.  In  the  first  part,  validation  of  indicial  theory 
is  performed  using  small  amplitude  steps  in  the 
boimdary  conditions  to  determine  the  indidal  responses 
of  the  flow  directly.  The  second  part  details  the 
implementation  of  indicial  theory  using  a  Laplace 
domain  approach  applicable  to  large  amplitude  (i.e., 
several  degrees)  motion,  including  cases  where 
amplitude  non-linearity  is  important.  In  the  third  part  of 
the  results,  the  accuracy  of  the  proposed  method  is 
evaluated  by  inverting  back  to  the  time  domain,  and  is 
applied  to  angular  motions  between  a  =  10*  and  a  =  15* 
at  a  Reynolds  number  Re  =  12,000. 


Part  1:  Direct  Validation,  Re  =  1,000. 

For  the  purpose  of  the  present  paper,  only  a 
limited  set  of  conditions  are  considered,  involving  small 
oscillations  about  two  angles  of  attack:  a  =  16*  and 
a  =  30*.  The  Mach  and  Reynolds  numbers  (respectively 
M„  =  0.2  and  Re  =  1,000)  are  fixed,  and  two  types  of 
indicial  responses  are  computed.  The  conventional 
determination  of  an  indicial  function  involves 
performing  a  step  in  the  boundary  condition  (in  this 
case,  a)  and  recording  the  response,  say  AOt  or  AO .  For 
t  >  0,  A0t(t)/Aa  is  the  indicial  response  of  ©t  with  respect 
to  its  boundary  condition  a  in  the  infinitesimal  limit 
Aa  0.  For  the  case  of  a  body  in  motion,  it  is  more 
relevant  to  consider  the  dynamic  indicial  response, 
calculated  as  the  Frechet  derivative®'^®  of  0,  with  respect 
to  a: 

0(o({)»t,T)  =  11m  (5) 

.  11a  .♦ttfl(O^H(g-T)Aa]  -  ♦Ja(^)l  j 

A  ^ 

Aa-»0 

where  the  step  in  boundary  condition,  Aa,  is  applied  at 
some  time  t  =  T,  and  H  designates  the  Heaviside  step 
function.  The  simple  case  of  a  ramp-and-hold  motion  is 
schematically  illustrated  in  Fig.  9.  For  example,  to  obtain 
the  dynamic  indicial  response  on  the  lift,  two  lift 
responses,  say  to  motions  (a)  and  (b)  in  Fig.  9,  are 
differenced  and  normalized  by  the  magnitude  of  the 
step,  Aa.  The  difference  between  maneuvers  (a)  and  (b) 
(at  16»)  or  (c)  and  (d)  (at  30*)  will  be  referred  hereafter 
as  the  "indicial  maneuver"  or  "indicial  motion."  In  this 
Part  1  of  the  validation  results,  the  step  magnitude  is 
small,  consistent  with  the  theory,  and  typically  on  the 
order  of  one  tenth  of  one  degree,  or  less. 

Hgure  10  compares  indicial  theoretical  predictions 
(lines)  and  direct  numerical  simulation  results  (symbols) 
for  both  phase  and  amplitude  of  the  fluctuating 
tangential  vorticity  flux.  The  direct  numerical 
simulations  involve  ramping  the  airfoil  from  a  =  0*  to 
a  =  16*  at  a  constant  pitch  rate,  da*/dt  =  0.6.  When  the 
airfoil  reaches  the  end  of  the  ramp,  it  starts  oscillating  at 
a  reduced  angular  frequency  u.  The  phase  (A#)  and 
amplitude  (A)  of  the  tangential  flux  response  at 
2%  chord,  A0t(t),  are  measured  after  5  or  6  cycles  of  the 
oscillation  (to  eliminate  transient  effects).  Both  A^  and  A 
are  defined  relative  to  their  quasi-static  coimterpart.  The 
lines  represent  indicial  theoretical  predictions,  assummg 
a  single  frequency  oscillatory  input,  £(t)  =  Aa  x  e'**^. 
The  dotted  line  corresponds  to  the  prediction,  Eq.  (3), 
using  the  dynamic  indicial  response.  The  solid  line 
corresponds  to  the  same  theoretical  prediction  using  the 
static  indicial  response  (dashed  line  in  Fig.  9).  The 
symbols  are  the  result  of  direct  numerical  simulations  for 
airfoil  harmonic  oscillations  about  a  =  16*.  The  circles 
are  associated  with  ramp-and-oscillate  maneuvers,  while 
the  triangle  corresponds  to  the  case  of  an  oscillation 
started  from  a  steady  angle  of  attack.  For  both  types  of 
maneuvers,  the  comparison  is  seen  to  be  excellent. 

A  similar  comparison  (a  =  16  *)  is  given  in  Fig.  11 
for  the  phase  and  amplitude  of  the  vorticity 
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accumulation,  O*,  between  direct  numerical  simulations 
and  the  indicial  theoretical  prediction.  The  direct 
numerical  simulations  require  each  on  the  order  of  a 
half-hour  to  an  hour  on  a  CRAY  Y-MP;  whereas,  the 
indicial  theoretical  prediction  takes  30  seconds  on  a 
VAX-780  for  100  frequencies.  As  in  Hg.  10,  theory  and 
numerical  experiment  agree  well.  But  when  the  same 
experiment  is  repeated  (Figure  12)  for  a  hold-angle  of 
30*  ((c)  and  (d)  in  Fig.  9),  the  agreement  with  indicial 
theory  is  poor.  The  non-linearity  of  the  vorticity 
accumulation  manifests  itself  in  the  form  of  a  low 
frequency  modulation  superimposed  on  the  basic  forcing 
frequency  (see  inset  of  Fig.  12).  A  major  difference 
between  o  =  16*  (Figs.  10  and  11)  and  a  =  30«  (Fig.  12)  is 
that  the  flow  is  asjrmptotically  tmsteady  for  the  latter. 
At  the  low  Reynolds  number  considered  here,  the 
solution  at  o  =  16»  is  steady.  It  is  important  to  realize, 
however,  that  for  larger  Re3molds  numbers,  for  example 
Re  =  12,000  and  beyond,  the  flow  is  both  separated  (to 
varying  degrees)  and  unsteady  at  all  angles  of  attack. 
As3rmptotic  (large  time)  unsteadiness  per  se  does  not 
preclude  the  use  of  indicial  theory.  However,  it  was 
found  in  the  case  of  Fig.  12  that  the  indicial  function  itself 
(calculated  from  the  relative  response  between 
maneuvers  (c)  and  (d))  does  not  converge  to  steady  state, 
a  situation  to  which  indicial  theory  does  not  apply. 


from  the  physics  of  the  cases  to  be  predicted.  Thus,  a 
straight-forward  time  domain  application  of  the  theory, 
as  stated  in  Eqs.  (3)  to  (5),  carmot  be  realized  unless  the 
indicial  motion  approximates  an  infinitesimal  step.  To 
remedy  this  problem,  we  now  consider  the  Laplace 
domain  version  of  the  indicial  formalism. 


The  Laplace  transform  of  Eq.  (3)  is  given  by 


f(s)  =  sl{s)c(s)  (6) 

whete7(s),  £(s),...  designate  the  Laplace  transforms  of  7 , 
£,  etc.,  respectively.  Keeping  in  mind  that  /  and  e 
represent  changes  with  respect  to  a  base  line  maneuver 
(denoted  "1”),  Eq.  (6)  can  be  inverted  as  follows 


7(t)  =  p^lt) 


e^{B)-e^{B) 


(7) 


where  L  designates  the  Laplace  operator.  In  (7), 
maneuver  "1"  is  characterized  by  the  boundary  condition 
schedule  £i(t)  and  by  the  absolute  flow  response  gjft); 
maneuver  "2"  is  characterized  by  which 

differs  from  {Ci(t),gi(t)}  after  time  t  =  0  (see  Figure  13). 
Naturally,  one  recovers  the  ori^nal  formulation  for 


e2(t)  =  ei(t)  +  AaH(t)  (8) 


To  circumvent  this  difficulty  and  extend  the  range 
of  validity  of  indicial  theory,  a  method  was  sought  to 
increase  the  amplitude  of  the  "indicial  maneuver.  The 
rationale  behind  this  approach  can  be  described  as 
follows.  Since  the  chaotic  dynamics  which  "contaminate" 
the  indicial  responses  are  an  inherent,  insuppressible, 
feature  of  the  flow  at  large  Reynolds  numbers  and  large 
angles  of  attack,  a  possible  approach  is  to  attempt 
instead  to  minimize  the  pseudo-random  or  chaotic 
component  of  the  response.  This  can  be  achieved  by 
maximizing  the  deterministic  portion  of  the  flow 
response  through  the  use  of  large  amplitude  indicial 
motions.  Of  course  there  will  be  an  upper  limit  to  this 
approach,  as  one  gets  into  the  non-linear  regime.  In 
effect,  the  new  approach  (explained  below.  Part  2)  relies 
on  "finite  amplitude"  (say  on  the  order  of  several 
degrees),  rather  than  infinitesimal,  indicial  maneuvers. 
These  are  used  to  infer  indicial  responses  in  which  the 
non-deterministic  component  is  kept  to  a  minimum. 

Part  2:  Laplace  Domain  Approach. 

There  are  several  difficulties  associated  with  the 
calculation  of  finite  step  responses  using  CFD.  Most 
often,  CFD  codes  will  not  run  for  these  conditions,  except 
for  the  smallest  of  step  sizes  and,  even  so,  using  very 
small  time  steps.  Even  if  a  stable  solution  can  be 
obtained,  its  accuracy  is  dubious  because  such  changes  in 
the  boundary  conditions  are  not  consistent  with  the 
temporal  and  spatial  orders  of  accuracy  of  the  numerical 
methods  being  used.^®  Furthermore,  the  flow  responses 
can  be  non-physical,  introducing  for  example  transient 
shocks  and  piston  effects  which  are  absent  in  the  real 
flow.  As  mentioned  previously,  it  stands  to  reason  that 
indicial  theory  cannot  be  used  when  the  physics 
contained  in  the  indicial  response(s)  of  the  system  differ 


with  Ae  =  const.  In  this  case:  =  Aa/s  and, 

obviously:  Pg(t)  =  A/i(t)/A£,  which  is  the  analog  of  Eq.  (5) 
when  (1)  is  exact.  In  the  general  case,  however,  note  that, 
if  Eq.  (7)  is  used,  62  and  need  not  differ  by  a  constant 
(see  Fig.  13),  nor  does  their  difference  need  to  be 
infinitesimaL  which  is  the  basis  of  the  present  method. 

The  validation  of  the  Laplace  domain-based 
method,  Eqs.  (3)  and  (7),  is  first  shown  for  low  angles  of 
attack  with  trailing  edge  flow  separation  (Re  =  12,000, 
M.  =  0.2).  To  illustrate  how  the  method  works,  consider 
Figs.  14  and  15.  Both  maneuvers  smoothly  ramp-up  the 
angle  of  attack  by  a  whole  degree,  from  o  =  4*  to  o  =  5*. 
Figure  14  depicts  the  history  of  the  tangential  vorticity 
flux  for  the  airfoil  at  rest  (curve  "1"),  for  the 
motion  (curve  "2"),  and  for  the  differential  response  ("2"- 
"1").  The  motion  is  the  same  fifth  order  polynomial  ramp 
as  in  Refs.  20  and  21,  chosen  primarily  for  its  smoothness 
(zero  first  and  second  derivatives  at  the  end  points).  The 
motions  considered  in  Figs.  14  and  15  differ  only  by  the 
overall  length  of  the  ramp,  or  "ramp  time,"  denoted  T. 
Figure  14  depicts  the  case  of  a  "fast"  maneuver  (T  =  0.1), 
and  Hg.  15  that  of  a  maneuver  taking  place  over  a  period 
40  times  longer  (T  =  4).  Note  that  the  difference  curve 
("2"."1'')  in  Fig.  14  is  not  the  indicial  response,  but  a  so- 
called  characteristic  response  from  which  the  indicial 
response  can  be  inferred,  using  (7).  The  goal  is  to 
evaluate  whether  the  slow  response  (T  =  4)  can  be 
predicted  from  the  application  of  indicial  theory,  when 
using  the  indicial  function  inferred  from  the  fast 
maneuver  (Hg.  14). 

Figure  16  compares  the  result  of  the  direct 
numerical  simulations  to  the  Laplace  domain  indicial 
prediction.  Two  pairs  of  ctuves  are  shown  in  the  figiue. 
In  each  pair,  one  curve  represents  the  flux  Laplace 
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transform  from  DNS  and  the  other  its  prediction  by  the 
indicia!  method.  For  T  =  4,  the  prediction  is  based,  as 
previously  mentioned,  on  an  indicia!  function  inferred 
from  the  T  =  0.1  nnaneuver.  The  two  curves  are  virtually 
indistinguishable  from  one  another  over  the  entire 
Laplace  spectrum,  with  the  exception  of  the  largest 
values  of  a  (corresponding  to  the  short  time  scales).  The 
agreement  is  equally  good  for  the  second  curve  pair 
shown  in  Fig.  16,  which  illustrates  the  prediction 
accuracy  for  the  fast  maneuver,  based  on  an  indicia! 
function  inferred  from  the  slow  one  (Fig.  15).  This  is,  of 
course,  to  be  expected  for  a  system  which  would  be 
perfectly  linear. 

In  Figures  14  through  16,  the  motions  that  were 
considered  were  all  of  the  same  amplitude:  6a  =  1  *.  In 
contrast,  if  one  attempts  to  use  the  indicial  function 
inferred  from  T  =  0.1  (Fig.  14)  to  predict  a  large 
amplitude  maneuver  from  a  =  4*  to  o  =  8*,  the  quality 
of  the  prediction  deteriorates  significantly.  This  is 
indicated  in  Fig.  17  by  the  spread  between  DNS  and 
prediction  for  s  <  1.  Note  that  all  four  curves  in  Fig.  17 
are  parallel  to  each  other  at  small  a  and  are,  in  fact, 
proportional  to  1/s.  Therefore,  the  source  of  the 
discrepancy  for  large  time  is  caused  by  amplitude  non¬ 
linearity  of  the  static  response.  In  effect,  the  linear 
assumption  that  the  increment  in  steady  vorticity  flux 
between  o  =  4*  and  a  =  8*  is  four  times  that  between 
a  =  4*  and  ct  =  5*  breaks  down,  even  at  these  low  angles 
of  attack.  This  is  to  be  expected  for  any  angle  change 
other  than  infinitesimal  since,  even  for  an  attached  flat 
plate  boundary  layer,  the  vorticity  flux  is  proportional 
to  U|.  In  other  words,  the  flow  response  must  be  non¬ 
linear.  Given  these  focts,  the  question  that  the  remainder 
of  this  paper  attempts  to  answer  is:  to  what  extent  can 
the  formalism  of  indicial  theory  be  used  to  predict  the 
flow?  And,  if  necessary,  how  must  the  theory  be 
amended  to  handle  the  cases  of  interest? 

In  Fig.  18,  the  accurate  prediction  of  the  large 
amplitude  case  shown  in  Fig.  17  is  achieved  by  using 
indicial  theory  based  on  an  equivalent  indicial  response 
which  is  inferred  from  a  large  amplitude  ramp 
maneuver,  from  o  =  4*  to  o  =  8*  (T  =  0.4).  This  indicial 
maneuver  differs  in  two  ways  from  that  considered 
previously.  1)  the  ramp  time  is  T  =  0.4  (instead  of  T  =  0.1) 
because  of  limitations  of  the  code  (maximum  da*/ dt), 
and  2)  the  amplitude  of  the  motion  (Aa  =  4*)  is  now  the 
same  as  the  motion  to  be  predicted.  Making  the 
amplitudes  the  same  will,  evidently,  take  care  of  the 
previous  discrepancy  due  to  quasi-static  non-linearity. 
What  the  comparison  of  Fig.  18  shows,  however,  is  that, 
once  quasi-static  amplitude  non-linearity  is  taken  into 
account,  the  unsteady  phase  lags  are  predicted 
accurately  using  the  formalism  of  indicial  theory.  In 
other  words,  the  unsteady  component  of  the  flow 
behaves  almost  linearly  if  the  indicial  function  is 
appropriately  stretched  in  amplitude.  Since  the 
comparisons  of  Figs.  16  and  18  are  only  along  the  real 
axis  of  the  Laplace  domain,  they  are  to  l«  regarded  only 
as  an  indication  that  indicial  theory,  used  in  the  above 
manner,  can  predict  the  vorticity  dynamics  of  the  flow. 
The  real  proof  of  concept  must  involve  all  of  the  Laplace 


domain  information  (i.e.,  inversion  back  to  the  time- 
domain)  and  must  apply  to  the  case  where  leading  edge 
separation  is  present.  This  validation  exercise  is  the 
object  of  Part  3. 

Part  3:  Validation  of  the  Laplace  Domain-Based 
Approach,  Re=  12,000. 

In  this  section,  the  implementation  of  indicial 
theory  as  described  in  Part  2  is  attempted  for  airfoil 
pitching  motions  ranging  from  a  =  10*  to  o  =  15*.  For 
all  of  these  cases,  the  Reynolds  number  is  Re  =  12,000, 
and  the  freestream  Mach  number  is  M.  =  0.2.  Figure  19 
illustrates  the  angular  motions  under  consideration. 
These  are  characterized  by  the  total  ramp  time,  T,  in  the 
case  of  ramp-up  or  ramp-down  motions,  or  by  the 
frequency,  F,  in  the  case  of  oscillations  (an  equivalency 
between  the  two  can  also  be  defined  by  considering  the 
first  half  of  the  oscillatory  cycle  to  correspond  to  a 
smooth  ramp,  and  vice-versa).  The  maneuver 
characteristics  are  summarized  in  Table  1,  along  with  the 
airfoil  motions  considered  in  Part  2. 


oC) 

T 

k=>rFc/U. 

(daVdt)^ 

~r 

to  5 

0.1 

15.70 

0.327 

4.0 

0.39 

0.008 

4 

to  8 

0.4 

3.93 

0.327 

4.0 

0.39 

0.033 

10 

to  15 

0.5 

3.14 

0.327 

0.7 

2.24 

0.234 

2.0 

0.78 

0.082 

2.5 

0.63 

0.065 

8.0 

0.20 

0.020 

Table  1.  Summary  of  Airfoil  Motion  Conditions 


In  the  examples  considered  here,  the  various 
maneuvers  are  all  referenced  to  the  fastest  ramp-up 
(T  =  0.5)  motion  in  Fig.  19.  This  means  that  the  raw 
vorticity  flux  or  vorticity  acounulation  responses  are  all 
processed  so  that  the  relative  flow  response  is  defined 
with  respect  to  the  case  T  =  03.  For  example,  the  relative 
vorticity  accumulation,  6o’,  for  T  =  2.5  is  defined  as: 

6^’(t)  =  •t.j  J  ■  *1-05 

and  so  on.  It  is  this  relative  flow  quantity  that  the 
indicial  theory  attempts  to  predict,  in  response  to  the 
relative  change  in  angle  of  attack,  6a,  defined  as: 

6a(t)  =  a7_25  -  07-05 

The  relative  angle  of  attack  variation  is  the  input  to  the 
indicial  prediction,  and  is  shown  in  Fig.  20  (solid  line). 
The  other  lines  in  Fig.  20  represent  the  relative  motion, 
6a,  for  the  various  ramp  times  and  frequencies  shown  in 
Fig.  19. 

Before  discussing  the  results  of  the  prediction 
method,  an  example  of  the  computed  flowfields  being 
analyzed  is  shown  in  Fig.  21.  Figme  21  depicts  contours 
of  the  primary  vorticity  at  two  time  instants  during  a 
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typical  pitch-up  maneuver  (T  =  0.5).  The  flowfield 
solutions  are  rotated  in  the  same  reference  frame  for 
comparison  purposes  (the  equations  of  motion  are 
solved  ill  an  inertial  frame  of  reference  in  which  the 
airfoil  is  rotated).  Figure  21  illustrates  a  t)rpical  result, 
obtained  using  the  computational  grid  shown  in  Fig.  1. 
The  vorticity  contours  at  t  =  0.4  (80%  into  the  maneuver) 
suggest  the  presence  of  an  essentially  attached  shear 
layer  (dashed  lines)  at  the  leading  edge.  The  layer 
separates  and  rolls  into  a  vortex  some  time  later  (solid 
lines),  as  is  illustrated  by  the  topological  changes  seen  in 
Fig.  21. 

The  changes  in  the  integrated  vorticity  fluxes  and 
vorticity  accumulation  at  the  leading  edge  during 
vortical  lift-off  and  roll-up  have  been  examined  in  detail. 
Typical  prediction  accuracy  is  shown  in  Figs.  22 
through  24. 

Figures  22, 23,  and  24  compare  the  DNS  results  for 
vorticity  accumulation  (solid  line)  to  the  indicial 
theoretical  prediction  (dashed  line).  The  latter  was 
transformed  to  the  time  domain  by  using  a  numerical 
inverse  Laplace  transform  technique  based  on  Laguerre 
functions."  The  prediction  is  obtained  using  the  method 
described  in  Part  2,  where  the  indicial  function  is 
inferred  in  the  Laplace  domain  from  an  indicial  motion 
(indicial  "pulse,"  in  the  present  case)  corresponding  to 
the  difference  between  the  two  fastest  ramps,  T  =  0.5  and 
T  =  0.7.  This  indicial  motion  is  shown  as  the  short 
negative  pulse  in  Fig.  20.  Figure  22  compares  the 
vorticity  accumulation  between  DNS  and  prediction  for 
a  ramp-up  motion  with  T  =  2.5.  The  prediction  is 
reasonably  accurate  at  short  times  and  even  appears  to 
capture  some  of  the  details  of  the  flow  response  up  to 
t  =  2.5. 

At  large  times,  there  are  intrinsic,  pseudo-random, 
fluctuations  due  to  vortex  shedding  and  other  chaotic 
phenomena  in  the  flow.  Thus,  the  relative  resfwnse  is 
inherently  unsteady,  even  long  after  the  motion  has 
ceased.  This  form  of  unsteadiness  is  uncorrelated  with 
the  motion  itself.  Recall  that  the  predictions  shown  here 
are  for  single  realizations.  The  use  of  phase-averaging 
on  multiple  realizations  would  presumably  average  out 
the  non-deterministic  portion  of  the  signal  and  increase 
the  accuracy  of  the  pr^ction.  Note  that  the  prediction 
itself  is  based  on  a  single  realization  of  ramp  maneuvers 
T  =  0.5  and  T  =  0.7.  Again,  realization  averaging  would 
be  desirable  to  increase  the  accuracy  of  the  prediction. 
The  latter  can  only  predict  the  deterministic  unsteadiness 
caused  by  changes  in  the  boundary  conditions,  6a.  In 
the  present  case,  realization  averaging  was  too  expensive 
(using  CFD-generated  results),  and  was  simulated 
instead  by  filtering  out  the  fluctuating  component  of  the 
indicial  response  for  t  >  3.2  (recall  that  the  indicial 
motion  is  completely  over  (6o  =  0)  by  t  =  0.7). 

In  spite  of  the  pseudo-random  component  of  the 
fluctuation  response,  the  accuracy  of  the  prediction 
during  the  initial  portion  of  the  motion  is  clearly 
illustrated,  both  in  the  ramp  cases.  Figs.  22  and  23,  and  in 
the  oscillatory  case.  Fig.  24.  The  slight  undershoot  and 


overshoot  of  the  prediction  for  t  s  0.7  are  believed  to  be 
related  to  the  accuracy  of  the  inverse  Laplace  transform 
procedure,  which  is  itself  an  approximation.  The 
consistency  with  which  even  small  scale  fluctuations  are 
initially  captured  by  the  indicial  prediction  in  all  three 
figures  (22, 23,  and  24)  suggests  that  these  fluctuations 
are  deterministic  in  nature  and,  therefore,  can  be 
predicted  using  indicial  theory,  as  described  in  Part  2. 
Even  at  large  times,  the  indicial  theoretical  prediction 
appears  to  capture  the  essential  features  of  the  flow 
response.  For  example.  Fig.  24,  the  phase,  amplitude, 
and  mean  of  the  oscillatory  vorticity  accumulation 
appear  to  be  predicted  reasonably  well. 

CONCLUDING  REMARKS 

Direct  numerical  simulations  have  been  used  to 
investigate  whether  indicial  theory  could  be  applied  to 
the  prediction  of  certain  key  aspects  (vorticity  dynamics 
characteristics)  associated  with  dynamic  stall. 

Preliminary  results  at  low  Reynolds  number 
(Re  s  12,000)  and  low  Mach  number  suggest  that 
integrated  vorticity  fluxes  and  vorticity  accumulation  at 
the  leading  edge  are  essentially  smooth  functions  of  their 
boundary  condition.  Based  on  this  observation,  the 
applicability  of  indicial  theory  was  first  evaluated 
directly  in  the  time  domain  using  step  responses  to 
calculate  the  indicial  functions  of  the  flow.  For  small 
motion  amplitudes  and  high  frequency,  unsteady 
vorticity  fluxes  and  accumulation  were  found  to  be 
predicted  well  by  indicial  theory. 

For  finite  amplitude  airfoil  motions,  large  angles  of 
attack,  and/or  larger  Reynolds  numbers,  the  straight 
forward  application  of  indicial  theory  gave  poor 
predictions.  To  circumvent  some  of  the  problems 
encountered  under  these  conditions,  a  Laplace  domain- 
based  method  using  finite  amplitude  indicial  motions 
was  implemented  and  tested  for  linear  and  non-linear 
airfoil  motions.  For  large  amplitude  (Aa  =  5*)  non-4inear 
maneuvers,  it  is  shown  that  the  linear  formalism  of 
indicial  theoiy  can  be  retained,  provided  that  the  indicial 
functions  are  inferred  in  the  Laplace  domain  from  a 
lefeience  indidal  motion  having  the  same  amplitude.  In 
other  words,  the  unsteady  component  of  the  flow 
behaves  almost  linearly  if  the  indicial  function  is 
appropriately  stretched  in  amplitude. 
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HGURES 


Hg.  1  Illustration  of  Computational  Grid  used  at 
Re  =  12,000;  (a)  Trailing  Edge  Region;  (b) 
Leading  Edge  (grey  area  indicates  the  control 
volume  used  for  the  calculation  of  integrated 
tangential  and  normal  fluxes  of  vorticity). 
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Fig.  2  Velocity  Vector  Plot  Illustrating  a  Typical 
Numerical  Solution  for  Re  =  1,000,  da'*’/ dt  =  0.6, 
M«  =  0.2. 


Fig.  5  Vortidty  Contours  in  Tangendal/Normal  Space 
and  Physical  Space  (Inset)  at  time  tU^c  =  0.075 
After  the  Application  of  the  No-Slip  Boundary 
Condition  (M«  =  0.1,  o  =  45®,  Re  =  100,000). 


Fig.  3  Streamlines  Illustrating  the  Absence  of  Flow 
Separation  for  the  Initial  Flowfield  Used  in  the 
Model  Problem  (M„  =  0.1,  t  =  0.0). 


Fig.  4  Vortidty  Contours  in  Tangential/JJormal  Space 
and  Physical  Space  (Inset)  at  time  tU^c  =  0.020 
After  the  Application  of  the  No-Slip  Boimdary 
Condition  (M„  =  0.1,  o  =  45 ",  Re  =  100,000). 
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Fig.  6  Comparison  of  Surface  Flow  Topology  at 
Re  =  100,000  (Top)  and  Re  =  1,000  (Bottom)  for 
a  =  30"  and  M„  =  0.1  (lines  demarcate  changes 
in  the  flow  direction  along  the  surface;  leading 
edge  at  s/s^ax  =  O-®;  arrows  indicate  flow 
direction). 
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Fig.  7  Comparison  of  Surface  Flow  Topology 
Diagrams  at  Re  =  1,000  Between  o  =  23», 
o  =  30*  and  a  =  45*  for  M„  =  0.1  (lines 
demarcate  changes  in  the  flow  direction  along 
the  surface;  leading  edge  at  s/s^^  =  0.5;  arrows 
indicate  flow  direction). 


Temporal  Evolution  of  the  Integrated 
Tangential  Vorticity  Flux  (Top),  Integrated 
Normal  Vorticity  Flux  (Center),  and  Flux  Angle 
(Bottom)  of  the  Leading  Edge  Vorticity 
(M.  =  0.1,  Re  =  1,000,  o  =  45-).  (Arrows 
indicate  direction  of  increasing  coordinate 
(streamwise  coordinate  in  top  graph;  normal 
coordinate  in  center  and  bottom  graphs). 
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Fig.  9  Schematic  Representation  Illustrating  the  Basic 
Procedure  Involved  in  Calculating  the  D)mamic 
Indidal  Fimction. 


DNS  Prediction 


A -  Static  Case 

O .  D)mamic  Case 


Fig.  11  Comparison  Between  Direct  Numerical 
Simulation  (Symbols)  and  Indicial  Theoretical 
Prediction  of  Phase  (Dashed  Line)  and 
Amplitude  (Solid  Line)  for  the  Accumulated 
Vortidty,  O*,  at  Re  =  1,000,  =  0.2,  a  =  16  • . 
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Fig.  12  Comparison  Between  Direct  Numerical 
Simulation  (Symbols)  and  Indicial  Theoretical 
Prediction  of  Phase  (Dashed  Line)  and 
Amplitude  (Solid  Line)  for  the  Accumulated 
Vortidty,  O*,  at  Re  =  1,000,  =  0.2,  a  =  30*. 

The  Inset  in  the  Bottom  Figure  Illustrates  the 
Presence  of  a  Low  Frequency  Modulation  of  the 
Output  at  a  =  30  •. 


Fig.  10  Comparison  Between  Indicial  Theoretical 
Prediction  (Lines)  and  Direct  Numerical 
Simulations  (Symbols)  for  the  Fluctuating 
Tangential  Vortidty  Flux  at  Re  =  1,000,  M«  =  0.2, 

a  =  16*.  ( - Static  Indicial  Response; 

•  .  •  •  Dynamic  Indicial  Response 
(da+/dt  =  0.6)). 
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MOTION  RESPONSE 


Fig.  13  Schematic  Illustrating  the  Use  of  Finite 
Amplitude  Indicial  Motions  to  Infer  the  Indicial 
Function  in  the  Laplace  Domain. 


0-4* ->5* 


Fig.  14  Illustration  of  the  Procedure  Used  to  Infer  the 
Indicial  Function:  Integrated  Tangential 
Vorticity  Flux  at  2%  chord.  Re  =  12,000,  T  =  0.1. 


0-4* ->5' 


Fig.  15  Illustration  of  the  Procedure  Used  to  Infer  the 
Indicial  Function:  Integrated  Tangential 
Vorticity  Flux  at  2%  chord.  Re  =  12,000,  T  =  4. 
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Fig.  16  Laplace  Domain  Comparison  of  Indicial 
Theoretical  Prediction  Versus  Direct  Numerical 
Prediction,  Re  =  12,000,  o  =  4  •  -»  5  • . 


Hg.  17  Laplace  Domain  Comparison  of  Indicial 
Theoretical  Prediction  Versus  Direct  Numerical 
Prediction  Illustrating  the  Deterioration  of 
Accuracy  for  Large  Time  in  the  Presence  of 
Quasi-Static  Non-Linearity,  Re  =  12,000, 
0  =  4’’  -♦8'. 
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Fig.  18  Laplace  Domain  Comparison  of  Indicial 
Theoretical  Prediction  Versus  Direct  Numerical 
Prediction  Illustrating  the  Recovered  Accuracy 
of  the  Prediction  When  Inferring  the  Indicial 
Function  from  the  Appropriate  Response, 
Re  =  12,000,0  =  4'' -»8«. 
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Fig.  19  Angle  of  Attack  Histories  Considered  for  the 
Validation  of  the  Laplace  Domain-Based 
Characteristic  Function  Approach  (T:  Ramp 
Time,  F:  Frequency). 


Fig.  22  Theoretical  Prediction  Versus  Direct  Numerical 
Simulation  for  Vorticity  Accumulation  During 
Ramp-Up  Motion,  10®  s  a  s  15*  (T  =  2.5), 

Re  =  12,000,  M„  =  0.2  ( - ARC2D, - 

Theoretical  Prediction). 


5a  (®) 


6a  =  a(t)  -  Oref^t). 


Fig.  21  Comparison  of  Primary  Vorticity  Contours 
Computed  at  Two  Instants  in  Time  During 
Ramp  Motion,  10®  s  a  s  15®  (T  =  0.5), 
Re  =  12,000,  M„  =  0.2. 


Fig.  23  Theoretical  Prediction  Versus  Direct  Numerical 
Simulation  for  Vorticity  Accumulation  During 
Ramp-Up  Motion,  10®  s  o  s  15®  (T  =  8.0), 

Re  =  12,000,  M„  =  0.2  ( - ARC2D, - 

Theoretical  Prediction). 

k= 0.785 


Fig.  24  Theoretical  Prediction  Versus  Direct  Numerical 
Simulation  for  Vorticity  Accumulation  During 
Oscillatory  Motion,  10®  s  a  s  15®  (F  =  0.25), 

Re  =  12,000,  M„  =  0.2  ( - ARC2D, - 

Theoretical  Prediction). 
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ABSTRACT 

It  has  been  suggested  that  eruptive  plumes  of 
vorticity  might  play  a  critical  role  in  the  physics  of 
vortex  formation  during  leading  edge  d3mamic  stall. 
Numerical  simidations  at  low  Reynolds  number  do  not 
seem  to  adequately  predict  this  phenomenon.  To 
explore  the  possibility  of  a  "bifurcation"  in  Reynolds 
number,  we  investigate  in  detaU  the  scaling  of  incipient 
laminar  separation,  vortex  formation,  and  shedding  with 
respect  to  Reynolds  number.  Numerical  simulations  are 
used  to  study  a  model  problem  in  which  a  two- 
dimensional  airfoil  remains  stationary  at  angle-of-attack, 
but  for  which  the  leading  edge  flow  separates  as  a  result 
of  an  impulsively  applied  no-slip  boundary  condition. 
The  calculations  are  laminar  (50,000  s  Re  s  400,(XK)),  and 
are  performed  for  a  =  15*  and  =  0.2.  The  resulting 
surface  flow  topology  is  analyzed  as  a  function  of 
Re5molds  number.  The  results  obtained  thus  far  appear 
to  contradict  the  h3rpothesis  that  a  form  of  bifvucation 
takes  place  at  some  intermediate  laminar  Reynolds 
number.  Furthermore,  times  and  locations  for  the  onset 
of  separation  bubbles,  vortex  formation,  and  feeding 
sheet  rupture  are  found  to  scale  according  to  various 
power  laws  of  the  Reynolds  number.  Re'",  with 
0.11  SOS  0.45. 


NOMENCLATURE 

Symbols  and  Abbreviations 
c  Airfoil  chord 
const.  Constant 

CFD  Computational  Fluid  Dynamics 
Cf  Skin  friction  coefficient  (Cf  =  21,^/ (P«Ul)) 

D  Dilatation  (D  =  du/  dx  +  dv/dy) 
dy  Normal  grid  spacing 

J  Jacobian(J=  au/dxav/ay- du/ay  av/ax) 

L.E.  Leading  Edge 

Freestream  Mach  number 
p  Pressure 

PS  Primary  separation  (same  as  PSD 
PSl  First  primary  separation 
PS2  Second  primary  separation 
R  Rupture 

Re  CaiordRe5molds  number  (Re  =  U.c/vJ 
s  Curvilinear  surface  coordinate 

551  First  secondary  separation 

552  Second  secondary  separation 

'Chief  Scientist,  Member  AIAA 
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553  Third  secondary  separation 

554  Foiuth  secondary  separation 

t  Time,  measured  from  the  instant  the  no-slip 

boundary  condition  is  applied 
TSl  First  tertiary  separation 
TS2  Second  tertiary  separation 
u  X-component  of  velocity 
U„  Freestream  velocity 

V  Y-component  of  velocity 
Vtx  Vortex 

X  Coordinate  aligned  with  the  airfoil 

y  Normal  coordinate 

a  Angle  of  attack 

&  Boundary  layer  thickness 

A  Variation 

Ax  Circumferential  grid  spacing 

V  Fluid  kinematic  viscosity 
p  Fluid  density 

T„  Wall  shear  stress 
n  Vorticity 

Subscripts 

e  Evaluated  at  the  outside  edge  of  the  boundary 
layer 

F  Vortex  formation 
inf  Freestream  value 
LE  Leading  edge 
max  Maximum 
min  Minimmn 
n  Normal 
o  Initial  condition 
t  Tangential 
w  Evaluated  at  the  wall 

n  Based  on  vorticity 
oo  Freestream  value 


INTRODUCTION 

In  recent  years,  research  into  model  problems  of 
unsteady  flow  separation^'®  has  shown  the  propensity  of 
boundary  layers  to  form  eruptive  plumes  of  vorticity. 
Their  connection  to  the  dynamic  stall®'^^  process  at  high 
Reynolds  numbers  is  discussed  in  Ref.  13.  Specifically, 
two  types  of  eruptions  are  distinguished.  The  first  one  is 
a  vortex-induced  eruption.  This  type  is  associated  with 
the  secondary  separation  provoked  by  the  presence  of  a 
primary  vortex  and  with  the  subsequent,  nonlinear, 
induction  which  results  from  this  interaction.  This 
eruption  process  has  been  observed and  computed 
in  a  number  of  flows.  The  second  type  of  eruption  is 
conjectured'®  to  be  a  spike-like  boundary  layer  response 
to  the  pressure  gradient.  This  spike  would  presumably 
occur  on  the  upstream  side  of  the  primary  separation 
bubble,  allowing  boundary  layer  vorticity  to  spill  away 
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from  the  surface  and  eventually  roll-up  into  the  primary 
vortex.  Unfortunately,  numerical  simulations  at  low 
Reynolds  number  do  not  seem  to  be  able  to  predict  this 
second  type  of  eruption.'^  Numerical  calculations  are 
possibly  too  diffusive,  or  cannot  resolve  some  of  the 
critical  sub-scale  physics.  **'^'*^  Another  possibility  is 
that  the  Re5molds  number  must  be  sufficiently  large  for  a 
sharp  eruption  to  take  place.*^  The  investigation  of  the 
existence  of  a  "bifurcation"  in  Reynolds  number,  i.e.,  the 
possibility  that  a  qualitative  change  in  the  behavior  of  the 
flow  might  occur  as  the  Reynolds  number  is  increased,  is 
the  subject  of  this  study. 

OBJECTIVE  AND  APPROACH 

The  objective  of  this  work  is  to  investigate  in  detail 
the  scaling  of  incipient  laminar  separation,  vortex 
formation,  and  shedding  with  respect  to  Reynolds 
number.  The  study  addresses  the  d)mamics  of  vorticity 
at  the  leading  edge  of  a  two-dimensional  airfoil  and  has 
application  to  dynamic  stall.  The  specific  goals  are  to 
analyze  the  above  scaling  properties  in  terms  of 
topological  changes  of  the  flowfield,  both  volumetric  and 
on  the  surface  of  the  airfoil. 

The  method  of  approach  is  to  use  numerical 
simulations  of  a  model  problem,  in  which  the  airfoil 
remains  stationary  at  angle-of-attack,  but  for  which  the 
leading  edge  flow  separates  as  a  result  of  impulsively 
applying  the  no-slip  boundary  condition.  This  model 
problem  is  well-suited  for  the  desired  investigation  of 
Reynolds  number  effects  because  the  initial  conditions 
are  independent  of  Reynolds  number  and  because  the 
temporal  and  spatial  scales  involved  are  small.  These 
characteristics  contribute  to  substantial  savings  in 
computer  resources,  as  compared  to  full  calculations  of 
dynamic  stall  at  "high"  Reynolds  number  (say. 
Re  1  100,000).  The  accuracy  and  reliability  of  current 
turbulence  models  is  insufficient  for  flow  physics 
investigations;  therefore,  the  present  investigation  is 
limited  to  "high"  Reynolds  number  laminar  simulations 
of  the  two-dimensional  Navier-Stokes  equations. 

No  experimental  data  is  available  for  this  model 
problem.  Therefore,  the  accuracy  of  two-dimensional 
Navier-Stokes  simulations  can  only  be  validated  by 
means  of  careful  grid  convergence  studies  and  by 
comparisons  between  different  simulations  of  the  same 
case.  Two  codes  are  used  for  this  comparison;  a  version 
of  NASA  Ames's  ARC2D  Navier-Stokes  solver*^  using  a 
time-implicit  integration  mode  and  second-order  spatial 
accuracy,  and  a  specially  developed  time-explicit  eighth- 
order  accurate  version  referred  to  as  ARC2D/8E. 


RESULTS 

Model  Problem 

The  study  of  a  model  problem  was  considered,  in 
which  the  airfoil  is  maintained  at  a  fixed  angle  of  attack 
and  the  no-slip  boundary  condition  is  impulsively 


applied.  The  surface  boundary  conditions  are  a 
prescribed  tangential  slip  velocity  distribution  u,(s)  for 
t  <  0,  and  zero  velocity  for  t  »  0.  The  tangential  slip 
velocity  at  t  <  0  is  determined  from  the  requirement  that 
=  0  everywhere  at  the  surface.  These  boundary 
conditions  result  in  attached  flow  at  t  =  0  which  has  the 
property  of  being  independent  of  Reynolds  number. 
Subsequent  to  the  application  of  the  no-slip  condition,  an 
unsteady  flowfield  develops.  This  flowfield  is 
characterized  by  the  complete  absence  of  trailing  edge 
separation  over  the  time  scales  of  interest,  and  the 
existence  of  a  well-defined  time  origin.  These  factors  are 
beneficial  in  terms  of  studying  the  initial  stages  of 
leading  edge  dynamic  stall,  including  scaling  effects  with 
Reynolds  munber. 

Based  on  our  previous  (impublished)  work  on  this 
model  problem,  at  least  two  distinct  stages  of  the 
incipient  separation  process  are  found.  Stage  I 
corresponds  to  the  formation  of  an  initially  symmetric 
(bubble-shaped)  separation  of  the  thin  primary  vorticity 
layer  which  lifts  away  from  the  surface.  This  stage 
appears  to  exhibit  qualitatively  similar  behavior  between 
Reynolds  numbers  over  the  investigated  range  (up  to 
5(X),0(X),  based  on  chord),  and  ends  with  the  beginning  of 
asymmetry  (lopsided  bubble),  first  evidenced  by  the 
accumulation  of  vorticity  in  the  reattachment  region  of 
the  separated  vorticity  layer. 

In  Stage  D,  the  accumulated  vorticity  develops  into 
a  vortex  which  subsequently  induces  a  sequence  of  rapid 
events,  the  first  of  which  is  the  formation  of  a  secondary 
recirculation  region  at  the  surface.  Stage  II,  in  the 
present  context,  is  defined  as  ending  with  the  "rupture" 
of  the  feeding  vortex  sheet,  and  the  subsequent 
convection  of  the  primary  dynamic  stall  vortex. 

Of  the  two  stages  distinguished  above.  Stage  II  is 
the  most  difficult  one  to  compute,  and  issues  of 
adequacy  of  the  numerical  resolution,  diffusion  of  the 
vortices  and  steepness  of  the  induced  pressure  gradients 
are  critical.  At  high  Reynolds  numbers,  a  second-order 
method  is  extremely  expensive  if  it  is  to  capture  the 
essential  physics.  The  special  feature  of  ARC2D/8E  is 
the  inclusion  of  a  high-order  accuracy  scheme  which 
permits  the  extension  of  the  model  problem  studies  to 
Reynolds  numbers  considerably  larger  than  previously 
possible  (Re  =  100,000  and  up). 

ARC2D/8E  code 

The  results  presented  here  are  obtained  using  the 
high-order  code  ARC2D/8E.  The  computational  method 
consists  of  central,  eighth-order,  finite-difference  spatial 
discretization,  and  an  explicit  third-order  Runge-Kutta 
time  integration  scheme.  Eighth-order  accuracy  is  also 
implemented  for  the  calculation  of  the  convective  fluxes, 
metric  terms,  and  artificial  dissipation.  Viscous. fluxes 
are  computed  with  fourth-order  accuracy.  The  order  of 
the  central  differencing  scheme  degenerates  to  second- 
order  at  the  wall.  The  high-order  artificial  dissipation 
algorithm  does  not  affect  modes  in  the  solution  which 
are  accurately  resolved  by  the  discretization  scheme. 
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The  accuracy  of  the  algorithm  has  been  validated  in 
various  studies.^®'^^  The  advantage  of  the  high-accuracy 
scheme  is  that  true  accuracy  can  be  attained  for  roughly 
half  of  the  grid-resolved  wave  number  range.  In 
comparison,  a  second-order  method  is  only  accurate 
over  approximately  5  to  10%  of  the  spectral  range,  while 
spectral  methods  typically  resolve  approximately  two 
thirds  of  the  theoretical  range.  Thus,  only  a  few  points 
(on  the  order  of  five)  are  required  to  resolve  a  shear  layer 
with  an  eighth-order  method.  Such  considerations  are 
critical  when  attempting  to  capture  the  physics  of  high 
Reynolds  ntimber  flows. 

Flow  Conditions  and  Grid  Convergence  Studies 

Four  Reynolds  munbers  have  been  investigated  in 
this  study;  these  are  (based  on  chord):  50,000, 100,000, 
200,000,  and  400,000.  All  simulations  are  for  a  fixed 
angle  of  attack  a  =  15*,  and  fixed  freestream  Mach 
number  M.  =  0.2.  The  results  shown  here  are  all 
computed  for  a  NACA0012  airfoil  with  a  rounded 
trailing  edge  (radius  of  curvature  r/c  =  0.0092).  This 
ensures  that  the  CFL  requirements  with  the  present  O- 
grid  are  not  dictated  by  the  trailing  edge  region.  The 
computational  grids  used  in  the  present  study  are  highly 
clustered  in  the  leading  edge  region  of  the  airfoil.  An 
example  of  a  computational  grid  used  in  this  study  is 
shown  in  Fig.  1.  For  this  grid,  476  points  were  used  in 
the  circumferential  direction,  340  of  which  were  placed 
in  the  first  4%  of  chord  on  the  upper  surface  of  the 
airfoil.  For  the  grid  convergence  studies,  the  grid  sizes 
ranged  from  224  to  476  in  the  circumferential  direction, 
and  between  60  and  83  for  the  normal  direction,  with  the 
outer  boundary  located  approximately  seven  chord 
lengths  away  from  the  airfoil. 

The  effect  of  circumferential  grid  spacing  is 
illustrated  in  Figs.  2  and  3.  Skin  friction  was  used  as  a 
sensitive  indicator  of  the  accuracy  of  the  solutions. 
Figure  2  (Re  =  200,000,  tU/c  =  0.084)  illustrates  the 
collapse  of  the  calculated  skin  friction  profiles  in  the 
region  of  interest  as  the  minimum  circumferential 
spacing,  Ax^^^/c,  is  reduced  from  l.OxlO"^  to  1.25x10*^. 
The  minimum  circumferential  grid  spacing  is 
representative  of  the  grid  in  the  critical  region 
0.01  s  x/c  s  0.04,  where  the  grid  was  designed  to  be 
nearly  imiform  in  the  streamwise  direction.  The  last  two 
grid  resolutions  (Ax^^^/c  =  2.5xl(H  and  1.25x10^)  yield 
virtually  identical  results.  The  value  AXjj,in/ c  =  1.25x10 
was  retained  for  Re  =  200,000  and  Re  =  400,000,  while 
AXn,in/c  =  2.5x10-*  and  Isx^iJc  =  5x10-*  were  used  for 
Re  =  100,000  and  Re  =  50,000,  respectively.  The  test  of 
Fig.  2  is  considered  to  be  fairly  stringent,  as  the 
resolution  requirements  for  the  collapse  of  the  Cf  profiles 
become  more  demanding  at  large  times.  This  form  of 
sensitivity  to  initial  conditions  is  illustrated  in  Fig.  3  at 
four  different  chordwise  positions.  Again,  it  is  seen  that 
grid  convergence  is  attained  for  Ax„,jn/ c  =  2.5x10-* . 

The  effect  of  normal  grid  spacing  at  Re  =  200,000  is 
shown  in  Fig.  4.  In  this  figure,  the  "error"  is  computed  as 
the  L-1  norm  of  the  difference  between  the  skin  friction 
distribution  and  the  "converged"  skin  friction 


distribution  at  some  fixed  time  instant  during  Stage  I. 
The  "converged"  Q  distribution  corresponds  to  the  best 
grid  resolution,  characterized  by  a  normal  grid  spacing 
at  the  wall  of  dy^/c  =  1.25x10-5.  Also  indicated  in  Fig.  4 
are  the  characteristic  slopes  associated  with  formal  third- 
and  fourth-order  accuracy.  The  data  presented  here  is 
seen  to  lie  somewhere  between  the  two,  which  is  not 
surprising  since  in  Stage  I  the  process  is  expected  to  be 
viscous-dominated,  and  since  fourth-order  accuracy  is 
the  maximum  achievable  for  the  viscous  fluxes  with  the 
present  code.  In  general,  lower  than  fourth-order 
accuracy  is  expected  near  the  wall,  since  the  formal 
accuracy  of  the  differencing  schemes  degrades  as  one 
approaches  the  boundary.  Based  on  the  result  of  Fig.  4,  a 
value  of  dyig/c  =  2.5x10-5  was  selected  for  Re  =  200,000 
and  400,000.  The  same  value  was  retained  for 
Re  =  100,000,  whfie  a  value  of  dy,e/c  =  5x10-5  was  used 
for  Re  =  50,000.  With  these  choices  of  parameters,  the 
normal  grid  spacing  at  the  wall  in  the  first  10%  of  chord 
is  better,  in  terms  of  viscous  units,  than  =  1  on 
average,  and  between  0  and  4  instantaneously. 

Anatomy  of  the  Flowfield 

A  typical  time  sequence  of  the  developing  vorticity 
field  in  shown  in  Fig.  5.  At  the  very  initial  stages,  a 
Stokes  layer  forms  after  the  imposition  of  the  no-slip 
boundary  condition.  This  is  followed  by  the 
development  of  a  thin  recirculation  bubble  along  the 
surface,  the  leading  and  trailing  edges  of  which  are 
marked  "PS"  (primary  separation).  The  bubble  is  initially 
symmetric,  but  eventually  becomes  lopsided  (see,  for 
example,  t  =  0.054  in  Fig.  5).  This  marks  the  end  of 
Stage  I.  Shortly  after,  a  vortex  forms  some  distance  away 
from  the  wall  in  the  downstream  region  of  the  primary 
recirculation  zone.  This  vortex,  in  turn,  induces  a 
secondary  separation  bubble.  This  first  instance  of 
secondary  separation  is  demarcated  by  the  symbols  SSI, 
visible  for  1 1  0.069  in  Fig.  5.  As  the  separating  primary 
shear  layer,  or  vorticity  feeding  sheet,  experiences 
instability  and  roll-up,  multiple  vortices  form.  As  in  the 
case  of  the  primary  vortex,  these  induce  new  secondary 
separation  zones.  For  instance,  the  leading  and  trailing 
edges  of  the  second  secondary  separation  are  indicated 
in  Fig.  5  by  the  symbols  SS2.  These  multiple  secondary 
recirculation  zones  eventually  strengthen  to  the  point 
where  they  too  become  vortices  which,  depending  on  the 
Reynolds  number,  may  induce  a  tertiary  separation,  and 
so  on.  In  Fig.  5,  counter-clockwise  vorticity  is  indicated 
by  dashed  lines.  From  the  last  two  graphs  of  Fig.  5,  we 
note  that  there  is  a  point  in  time  at  which  mutual 
induction  between  the  primary  dynamic  stall  vortex  and 
the  first  secondary  (counter-rotating)  vortex  leads  to  an 
eruption  of  vorticity  away  from  the  wall.  This  eruption 
appears  to  serve  as  a  trigger  for  the  rupturing  of  the 
vorticity  feeding  sheet,  in  agreement  with  Refs.  16 
and  22.  After  the  vortex  sheet  rupture,  the  primary 
dynamic  stall  vortex  is  free  to  convect  away  from  the 
leading  edge,  as  may  be  seen  from  the  last  graph  of 
Fig.  5. 

A  comparison  of  the  pressure  fields  at  the  various 
Reynolds  numbers  is  given  in  Fig.  6  at  tU/c  =  0.102. 


3- 


Although  the  different  flowfields  are  not  at  the  same 
stages  of  development  (the  larger  the  Reynolds  number, 
the  more  advanced  the  solution),  the  presence  of  vortex 
cores  can  be  clearly  seen  from  the  existence  of  closed  low 
pressure  contours.  To  accurately  resolve  the  secondary 
flows  induced  by  these  complex  pressure  fields,  it  is 
important  to  capture  as  high  a  portion  of  the 
wavenumber  spectrum  as  possible.  Thus,  use  of  the 
high-order  accuracy  code  ARC2D/8E  is  essential  at  the 
higher  Reynolds  numbers.  As  previously  mentioned, 
ARC2D/8E  uses  eighth-order  linear  and  nonlinear 
artificial  dissipation.  At  the  relatively  high  Reynolds 
numbers  considered  in  this  study,  a  concern  is  the  level 
of  artificial  dissipation,  compared  to  the  "natural" 
viscous  dissipation  of  the  flow.  This  comparison  was 
performed  in  a  systematic  manner  at  all  Reynolds 
numbers.  For  all  cases  investigated,  the  total  flux 
contribution  due  to  artificial  dissipation  was  less  than 
0.1%  of  the  viscous  contribution,  in  all  vortical  regions  of 
the  flow.  This  ratio  progressively  increases  away  from 
the  body,  as  viscous  fluxes  become  negligible  and  the 
artificial  dissipation  increases  due  to  the  coarser  grid 
resolution. 

Figure  7  depicts  the  space-time  evolution  of  the 
skin  friction  in  the  leading  ^ge  region  of  the  flow,  for 
Re  =  50,000  and  Re  =  200,000.  Such  surface 
representations  are  found  to  correlate  well  with  the 
instantaneous  topological  structure  of  the  flow.  A 
noticeable  feature  is  the  presence,  particularly  at  the 
higher  Reynolds  numbers,  of  a  "latency"  period  (e.g.. 
Fig.  7, 0.078  s  tU/c  s  0.090,  Re  =  200,000),  during  which 
the  vortex  system  grows  in  strength,  with  little  change  in 
either  size  or  position.  This  is  a  result  of  a  nearly 
stationary  solution  of  the  equations  of  motion.  The  end 
of  this  latency  f)eriod  is  associated  with  the  "rupture"  of 
the  vorticity  feeding  sheet  and  the  release  of  the  primary 
vortex,  as  indicated  in  Fig.  7  by  the  change  in  slope  of  the 
skin  friction  contours.  Again,  the  phenomenon  is  more 
clearly  marked  at  the  higher  Reynolds  numbers, 
although  the  basic  phenomenon  is  believed  to  be  still 
present  at  Re  =  50,0(K).  At  the  lower  Reynolds  numbers, 
viscous  diffusion  tends  to  blur  some  of  these  effects. 
Most  of  the  calculations  performed  in  this  study  are 
stopped  shortly  after  the  end  of  Stage  II,  since  our  goal  is 
to  study  the  physics  of  incipient  vortex  formation  and 
dynamic  stall.  At  the  final  time  instant  of  the 
calculations,  the  released  vortex  has  traveled  a  short 
distance  from  the  remainder  of  the  attached  vortex 
system. 

The  space-time  skin  friction  diagram  of  Fig.  6  can 
be  simplified  somewhat  by  considering  only  the  zero 
contours  of  skin  friction.  In  this  manner,  the  resulting 
lines  delimit  adjacent  regions  of  forward  and  reverse 
flow.  This  permits  the  identification  of  primary, 
secondary  and,  if  applicable,  tertiary  separation  regions. 
The  resulting  surface  flow  topology  at  Re  =  100,000  is 
shown  in  Fig.  8,  where  the  symbols  PS,  SSI,  SS2,  TSl, 
and  TS2  designate  respectively  the  onset  of  primary 
separation,  "first"  secondary  separation,  "second" 
secondary  separation,  "first"  tertiary  separation,  and 
"second"  tertiary  separation,  as  identified  in  Fig.  5.  In 


addition,  "R"  indicates  the  time  of  vortex  sheet  rupture, 
corresponding  to  the  end  of  the  previously  mentioned 
latency  {)eriod  and  the  beginning  of  free  convection  of 
the  primary  dynamic  stall  vortex.  The  solid  s)rmbols 
indicate  the  trajectory  of  the  topological  "center" 
associated  with  the  instantaneous  streamlines. 
Reference  17  illustrates  how  critical  point  theory^  can  be 
used  as  a  sensitive  means  to  analyze  the  tojx)logy  of  the 
velocity  field.  In  particular,  this  work*'  was  able  to 
pinpoint  the  origin  of  the  primary  stall  vortex.  Using 
similar  concepts,  it  is  found  in  the  present  flowfields  that 
the  critical  points  are  "centers,"  characterized  by 
essentially  zero  dilatation,  D,  and  positive  Jacobian,  J. 
Immediately  after  the  onset  of  surface  flow  reversal  in 
Fig.  8  (marked  "PS"),  the  value  of  the  Jacobian  nears  zero 
and,  hence,  the  center  is  degenerate  and  its  topology  is 
associated  with  pure  shear.  Although  the  center  does 
not  become  a  vortex  per  se  until  approximately 
tU/c  =  0.057  (shortly  before  the  onset  of  secondary  flow 
separation  "SSI"),  it  is  noteworthy  that  its  trajectory 
correlates  well  with  the  downstream  edge  of  the  primary 
and  induced  secondary  bubbles. 

The  relation  between  the  topological  center  and  a 
vortex  was  elucidated  by  Perry  and  Chong.^^  In  the 
present  study,  the  point  at  which  the  center  is  associated 
with  the  core  of  a  vortex  was  found  to  be  most  clearly 
defined  by  the  pressiue  field.  For  a  fairly  long  portion  of 
time  during  Stage  I,  (tU/ c  s  0.05  in  Fig.  8),  the  vorticity 
field  topology  corresponds  to  a  lifted  vorticity  layer 
beneath  which  a  thin  region  of  reverse  flow  exists. 
During  this  potion,  the  location  of  the  topological  center 
"rides"  the  vorticity  line,  gradually  moving  from  front  to 
back  and  settling  at  about  two  thirds  of  the  bubble 
extent.  In  contrast  to  the  rather  benign  behavior  of  the 
vorticity  during  that  portion,  the  end  of  Stage  I  is 
marked  by  a  sudden  change  from  boundary  layer-like 
behavior  (negligible  normal  pressure  gradient)  to  vortex¬ 
like  behavior,  characterized  by  the  presence  of  closed 
pressure  contours  associated  with  a  low  pressure  core. 
At  later  stages  of  the  evolution  of  the  flow,  local  pressure 
and  vorticity  fields  act  to  reinforce  each  other.  The 
vorticity  accumulates  further  in  the  low  pressure  region, 
and  increased  vorticity  results  in  a  further  lowering  of 
the  pressure.  This  interaction  rapidly  results  in  a  formed 
vortex  and  marks  the  beginning  of  what  was  referred  to 
earlier  as  Stage  H.  At  that  point,  vorticity  and  pressure 
fields  associated  with  vortices  are  well-correlated.  An 
example  at  the  end  of  Stage  11  is  given  in  Fig.  9.  The  top 
graph  in  Fig.  9  illustrates  the  strong  induced  normal 
entrainment  of  vorticity  immediately  "behind"  the 
(second)  dynamic  stall  vortex  at  Re  =  400,000.  The 
bottom  graph  shows  the  corresponding  instantaneous 
pressure  field. 

The  connection  between  flowfield  topology  and 
vortex  formation  is  further  explored  in  Fig.  10.  In  this 
figure,  the  initial  temporal  increase  of  the  Jabobian  value 
at  the  critical  point  is  characteristic  of  the  emergence  of  a 
center,  from  pure  shear  to  "pancake"-like  to  nearly 
circular.  After  the  point  of  vortex  formation  (associated 
with  the  formation  of  a  low  pressure  core),  the  Jacobian 
increases  further,  with  marked  stages  coinciding  with 


4- 


the  onsets  of  secondary  and  tertiary  separation.  Most 
importantly,  the  evolution  of  the  formed  vortex  is 
characterized  by  a  distinctive  drop  in  core  pressure.  This 
feature  was  observed  at  all  of  the  investigated  Reynolds 
munbers. 

Reynolds  Number  Scaling 

A  comparison  of  zero  skin  friction  diagrams  for 
Re  =  100,000,  Re  =  200,000,  and  Re  =  400,000  is  given  in 
Fig.  11.  As  in  Fig.  8,  the  symbols  indicate  the  times  and 
locations  associated  with  the  onset  of  flow  separation 
regions,  and  are  generically  indicated  for  Re  =  400,000. 
Although  not  repeated  in  the  figiue,  the  same  symbols 
also  apply  to  the  topological  features  associated  with  the 
other  Re5molds  numbers.  Also,  note  the  presence  of  a 
second  primary  separation  (PS2).  The  latter  was  also 
present  at  Re  =  100,000,  but  is  not  visible  on  the  scale  of 
Figs.  8  and  11.  This  second  primary  separation  is 
associated  with  a  re-separation  of  the  boundary  layer 
downstream  of  reattachment.  As  the  Re)molds  number 
increases,  the  overall  scale  of  the  leading  edge 
bubble /vortex  system  diminishes,  and  the  reattached 
flow  still  experiences  an  adverse  pressure  gradient 
which  causes  it  to  undergo  a  mild  separation.  The 
evidence  gathered  thus  far  suggests  that,  over  the 
investigated  range,  this  second  primary  separation  is 
completely  independent  and  does  not  interfere  with  the 
leading  edge  vortex  system.  For  the  purpose  of  clarity. 
Re  =  50,000  is  not  shown  in  Fig.  11.  The  important  point 
illustrated  by  this  figure  is  the  remarkable  degree  of 
similarity,  including  details,  between  results  at  the 
Reynolds  numbers  studied.  This  similarity  is  further 
exploited,  below,  to  extract  basic  scaling  relationships. 

Instantaneous  vorticity  fields  for  Re  =  50,000, 
100,000, 200,000,  and  400,000  are  compared  in  Fig.  12. 
The  value  of  the  vorticity  contours  is  kept  identical 
between  the  four  flowfields.  The  time  instant  chosen  for 
this  representation  corresponds  approximately  to  the 
end  of  Stage  II,  i.e.,  the  "rupture"  of  the  feeding  vortex 
sheet  and  the  subsequent  release  of  the  dynamic  stall 
vortex.  The  vortex  systems  in  Fig.  12  differ  primarily  in 
strength  and  scale  but  exhibit  topological  similarity. 
This  figure  shows  the  relative  physical  scale  of  the 
vortical  flowfields  as  a  function  of  Reynolds  number. 
Their  qualitative  similarity  is  seen  more  clearly  in  the 
expanded  views  provided  in  Fig.  13. 

From  zero  skin  friction  plots  such  as  Fig.  11,  one 
can  determine  with  reasonable  accuracy  the  time  of  onset 
of  the  various  separation  bubbles,  as  well  as  their 
corresponding  location.  The  scaling  of  the  separation 
time  of  onset  with  respect  to  Reynolds  number  is 
summarized  in  Fig.  14.  The  lines  in  the  plot  correspond 
to  straightforward  power  fits,  based  on  the  four  data 
points.  The  time  of  incipient  separation  for  the  primary 
bubble  is  fovmd  to  scale  as  Re-"  “,  while  the  first,  second, 
and  third  occurrences  of  secondary  separation  are  found 
to  scale  as  Re-®-^®,  Re-®-^®,  and  Re"®-^,  respectively.  For 
reference,  the  Reynolds  number  scaling  of  the  upstream 
boundary  layer  is  indicated  in  Fig.  15.  It  is  clear  that  the 
boundary  layer  thickness  6q,  at  the  leading  edge  of  the 


airfoil,  scales  as  Re'^^^  within  measurement  uncertainty. 
Therefore,  the  various  separation  events  do  not  scale  on 
the  upstream  boundary  layer  thickness  alone,  as 
expected  from  previous  model  flow  separation  studies.^*® 

As  a  final  note,  in  addition  to  the  scaling  of 
incipient  bubble  separation,  the  times  and  locations  of 
vortex  sheet  "rupture,"  as  well  as  those  associated  with 
vortex  formation,  have  been  recorded  for  the 
investigated  Reynolds  numbers.  "Vortex  formation  is 
defined,  rather  imprecisely  at  this  time,  as  the  transition 
between  an  elongated  separation  bubble  (characterized 
primarily  by  counter-flowing  streams  separated  by  a 
shear  layer)  and  a  vortex  structure  (characterized  by 
significant  normal  momentum  transfer  and  the  existence 
of  a  low  pressure  core).  The  time  instants  corresponding 
to  the  primary  separation  onset  (PSl),  vortex  formation, 
and  shear  layer  rupture  are  also  shown  (Fig.  16)  to  scale 
approximately  as  simple  power  laws  of  the  Re)molds 
munber.  Specifically,  the  time  instant  corresponding  to 
vortex  formation  scales  as  Re"®"^^,  while  that  of  rupture 
approximately  scales  as  Re"®‘^®.  Similar  power  law 
scalings  are  found  to  govern  the  location  of  these  events, 
as  shown  in  Fig.  17. 

If  some  form  of  boundary  layer  eruption  leading  to 
the  formation  of  the  primary  d)mamic  stall  vortex  were  to 
take  place,  modified  scaling  properties  would  be 
expected  from  the  emergence  of  new  physics.  In 
particular,  the  vorticity  spill-over  associated  with  a 
tightly  focused  eruption  would  presumably  result  in 
early  vortex  formation  and  rupture.  Instead,  the  present 
results  indicate  a  smooth  scaling  behavior  and  the  well- 
defined  emergence,  from  the  lifting  vorticity  layer,  of  a 
topological  center  which  precedes  vortex  formation. 


CONCLUDING  REMARKS 

Preliminary  results  into  the  Reynolds  number 
scaling  of  two-dimensional  incipient  flow  separation, 
vortex  formation,  and  rupture  from  its  feeding  sheet 
indicate  the  existence  of  qualitatively  similar  behavior 
between  the  chord  Reynolds  numbers  of  50,000  and 
400,000.  The  results  obtained  thus  far  appear  to 
contradict  the  hypothesis  that  a  form  of  bifurcation  may 
take  place  at  some  intermediate  laminar  Reynolds 
number:  while  the  flow  solutions  obtained  at  Re  =  50,000 
are  naturally  more  dissipative  than  for  the  highest 
Reynolds  numbers  computed  in  this  study,  it  remains 
that  the  sequence  of  phenomena  leading  to  the  release  of 
a  simulated  dynamic  stall  vortex  bear  strong  qualitative 
similarities.  In  particular,  times  and  locations  of  onsets 
of  specific  dynamical  events  are  found  to  scale  according 
to  simple  power  laws  of  the  Reynolds  number.  Re'**,  with 
0.11  s  /j  s  0.45. 
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T)^ical  Temporcd  Evolution  of  the  Vorticity  Field  (Re  =  100,000) 
Arrows  Indicate  Locations  of  Skin  Friction  Sign  Reversal, 
(continued  from  previous  page) 


Fig.  8  Temporal  Evolution  of  Surface  Flow  Topology 

and  Primary  "Vortex"  Trajectory  at 
Re  =  100,000.  Lines  are  Zero  Skin  Friction 
Contours.  "R"  Denotes  the  Time  of  Rupture  of 
the  Vortidty  Feeding  Sheet 
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Rg.  10  Temporal  Evolution  of  Key  Features  Used  to 
Describe  the  Primary  Vortex  (O  Jacobian, 
■  Pressure  at  the  Core). 
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Fig.  11  Similarity  of  Surface  How  Topologies  across 
Reynolds  Numbers. 
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Fig.  12  Comparison  of  Vorticity  Fields  at  the  End  of 
Stage  II,  for  Re  =  50,000,  Re  =  100,000, 
Re  =  200,000,  and  Re  =  400,000.  (Scale  is 
Maintained  the  Same  Between  Plots). 
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Fig.  13  Expanded  View  of  the  Vorticity  Fields  Shown 
in  Fig.  12  Illustrating  Their  Qualitative 
Similarity  at  the  End  of  Stage  H. 
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Fig.  14  Reynolds  Number  Scaling  of  Time  of 
Appearance  of  Incipient  Primary  (PSD  and 
Secondary  (SSI,  SS2,  SS3,  SS4)  Separation 
Bubbles  (see  Hg.  11). 


He.  16  Rejmolds  Number  Scaling  of  Time  of  Vortex 
Formation  and  Time  of  Vortex  Sheet  Rupture. 
Corresponding  Spadal  Locations  are  Shown  in 
Fig.  17. 


Fig.  15  Inverse  Square  Root  Scaling  of  the  Upstream 
Boundary  Layer  Thickness  with  Reynolds 
Number. 


Fig.  17  Reynolds  Number  Scaling  of  Location  of 
Vortex  Formation  and  Location  of  Vortex  Sheet 
Rupture.  Corresponding  Temporal  Scaling  is 
Given  in  Fig.  16. 
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COMPLEMENT  TO  AIAA  PAPER  NO.  94-2339 


Several  results  not  shown  in  AIAA  Paper  No.  94-2339  are  summarized  here  to  clarify 
the  methods  used  and  provide  further  documentation  of  the  cases  computed  in  that  paper. 

Comparison  between  ARC2D/8E  and  ARC2D,  As  mentioned  in  Ref.  Cl,  p.  2,  no 
experimental  data  is  available  to  compare  against  the  model  semi-impulsive  problem. 
Therefore  the  accuracy  of  the  calculations  was  verified  by  two  means.  The  first  one  is 
through  careful  grid  convergence  studies.  The  second  was  to  compare  the  computed  results 
from  two  different  codes.  The  two  codes  that  were  used  in  this  study  are  NASA  ^^les  s 
ARC2D  Navier-Stokes  finite  difference  solver,  and  a  specially  developed  explicit  eighth- 
order  accurate  version  named  ARC2D/8E. 

A  summary  of  the  grid  convergence  studies  performed  at  Re  =  200,000  is  given  in 
Appendix  B.  These  studies  were  carried  out  using  the  ARC2D/8E  code.  Appropriate 
computational  meshes  for  Re  =  100,000  and  Re  =  400,00  were  conservatively  designed  from 
these  results  and  from  simple  scaling  laws  relating  grid  spacing  to  Reynolds  ^mber.  A 
separate  grid  resolution  study  for  the  Re  =  800,000  case  is  given  in  Appendix  D.  The  present 
section  concerns  an  exercise  which  was  not  reported  in  Ref.  Cl,  namely  the  comparison 
between  the  results  of  the  implicit  second-order  code  ARC2D  and  the  explicit  eighth- 

order  code  ARC2D/ 8E. 

Figures  C.l  through  C.3  compare  the  vortical  flowfield  computed  with  ARC2D  and 
ARC2D/8E  at  successive  times,  tU/c  =  0.030, 0.051,  and  0.060,  for  the  chord  Reynolds 
number  Re  =  400,000.  The  results  are  in  good  agreement,  both  qualitative  and  quantitative, 
at  short  times.  However,  the  ARC2D  solution  breaks  down  (Fig.  C.3)  at  later  times.  This 
behavior  is  consistent  with  the  appearance  of  fine  scales  in  the  flowfield,  the  (numencal) 
growth  of  which  cannot  be  contained  by  the  artificial  dissipation  scheme  used  m  ARC2D.  m 
contrast,  the  solution  remained  stable  throughout  the  calculation  with  ARC2D/8E.  The 
computational  method  used  in  ARC2D/8E  consists  of  central,  eighth-order,  finite-difference 
spatial  discretization,  and  an  explicit  third-order  Runge-Kutta  time  integration  scheme. 
Eighth-order  accuracy  is  implemented  for  the  calculation  of  the  convective  fluxes,  eighth- 
order  for  the  metric  terms,  and  eighth-order  for  the  artificial  dissipation.  Viscous  fluxes  are 
computed  with  fourth-order  accuracy.  The  order  of  the  central  difference  scheme 
degenerates  to  second  order  at  the  wall.  The  high-order  artificial  dissipation  algorithm  does 
not  affect  modes  in  the  solution  which  are  accurately  resolved  by  the  di^retization  scheme. 
The  accuracy  of  the  algorithm  has  been  validated  in  various  studies.  '  The  advantage  of 
the  high-accuracy  scheme  is  that  true  accuracy  can  be  attained  for  roughly  half  of  the  grid- 
resolved  wave  number  range.  In  comparison,  a  second-order  method  is  only  accurate  oyer 
approximately  5  to  10%  of  the  spectral  range,  while  spectral  methods  typically  resolve 
approximately  two  thirds  of  the  theoretical  range.  Such  considerations  are  critical  when 
attempting  to  capture  the  physics  of  high  Reynolds  number  flows,  and  are  believed  to  relate 
to  the  failure  of  ARC2D  at  large  times. 

From  the  results  shown  in  Figs.  C.l  and  C.2,  it  is  tempting  to  conclude  that  the  use  of 
ARC2D  might  be  justified  provided  that  one  is  only  interested  in  topological  changes  of  the 
flow  at  early  times.  Indeed,  if  it  can  be  argued  that  ARC2D  is  as  accurate  as  ARC2D/8E  for  a 
specific  purpose,  then  the  implicit  algorithm  (resulting  in  higher  Courant  numbers  and,  thus, 
less  computational  work)  would  seem  more  attractive.  In  reality  however,  ARC2D,  when 
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run  with  a  CFL  number  of  15  produced  poor  results  in  comparison  to  ARC2D/8E,  even  for 
the  prediction  of  the  time  of  first  flow  reversal  at  the  surface.  The  results  for  ARC2D  were 
found  only  to  begin  to  achieve  time-step  independence  for  CFL  numbers  lower  than  three. 
This  is  illustrated  in  the  surface  flow  topology  diagram  comparison  of  Fig.  C.4.  For 
reference,  ARC2D/8E  was  run  with  a  CFL  number  of  1.5.  Taking  into  account  the  extra 
costs  associated  with  matrix  inversion,  the  computational  savings  associated  with  running 
ARC2D,  rather  than  ARC2D/8E,  were  minimal  for  this  problem.  Therefore  ARC2D/8E  was 
the  program  of  choice  for  this  study.  In  Appendix  D,  results  computed  with  ARC2D/8E  for 
a  chord  Reynolds  number  Re  =  800,000  are  presented. 

From  the  ARC2D  calculations  with  (At)a/c  =  9.4x10'^  (solid  line  in  Fig.  C.4)  we  show, 
for  completeness,  the  computed  characteristic  times  associated  with  critical  changes  of  the 
surface  flow  topology  during  the  early  stages  of  flow  development.  These  are  indicated  as 
shaded  symbols  in  the  Reynolds  number  scaling  plot  of  Fig.  C.5.  The  open  symbols 
correspond  to  the  ARC2D/8E  runs  at  various  Reynolds  numbers,  from  Re  =  50,000  to 
Re  =  400,000.  As  in  Appendix  B,  the  notation  used  here  is  as  follows:  the  time  of  first 
occurrence  of  the  primary  flow  separation  is  denoted  "PSl."  Successive  times  for  the 
occurrence  of  secondary  flow  separation  regions  are  denoted  (in  their  order  of  appearance) 
"SSI,"  "SS2,"  "SS3,"  etc.  Similarly,  tertiary  separation  is  denoted  ’TSl,"  ’TS2,"  and  so  on.  The 
results  presented  in  Fig.  C.5  indicate  that  the  ARC2D  determinations  of  PSl,  SSI  and  SS2  fall 
within  the  uncertainty  range  of  the  power-law  fits  discussed  in  Appendix  B. 

Reynolds  Number  Scalinp.  Similarly  to  Fig.  C.5,  simple  scaling  relationships  were  given  in 
Appendix  B  for  the  times  of  primary  vortex  formation  and  primary  vortex  release  (or 
"rupture"  of  the  feeding  vorticity  sheet).  Figures  12  and  13,  p.  12,  of  Appendix  B  illustrated 
the  qualitative  similarity  of  the  vortical  flowfields  at  the  end  of  Stage  II.  This  qualitative 
similarity  is  carried  one  step  further  in  Fig.  C.6.  In  this  figure,  instantaneous  vorticity  fields 
for  Re  =  50,000,  Re  =  100,000,  Re  =  200,000,  and  Re  =  400,000  are  compared;  the  value  of  the 
vorticity  contours  is  kept  identical  between  all  four  flowfields.  The  time  instant  is  chosen  to 
be  on  the  order  of  15%  (±  2%)  beyond  the  point  of  rupture  of  the  feeding  vorticity  sheet.  A 
vortex-induced  eruptive  plume  of  vorticity  is  visible  at  the  three  higher  Reynolds  numbers 
depicted  in  Fig.  C.6.  In  all  likelihood  the  flowfield  at  Re  =  50,000  would  also  exhibit  the 
same  type  of  eruption  at  a  slightly  later  time.  This  is  based  on  the  overall  similarity  of  the 
flowfields  and  observations  at  higher  Reynolds  numbers  at  earlier  times. 

Further  details  of  the  temporal  evolution  of  the  flowfield  at  Re  =  400,000  (Figs  C.7 
through  C.15)  reveal  the  existence  of  multiple  vortex-induced  eruptions  (see,  e.g..  Fig.  C.9). 
After  the  end  of  Stage  II  (i.e.,  the  rupture  of  the  primary  vortex  feeding  sheet),  Reynolds 
number  similarity,  in  the  sense  described  in  Appendix  B,  is  not  observed.  In  addition  to 
multiple  ejections,  vortex  pairing  events  take  place.  Some  of  these  complex  vortical 
interactions  take  place  earlier  still  at  higher  Reynolds  numbers.  This  is  described  further  in 
Appendix  D,  for  the  Re  =  800,000  case. 
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Fig  C.2(a)  Comparison  of  Instantaneous  Vorticity  Contours  Computed  Using  ARC2D  (Top) 
and  ARC2D/8E  (Bottom),  tU/c  =  0,051,  Re  =  400,000.  Spacing  between  contours 
is  An  =  80;  dashed  lines  indicate  negative  levels).  Continued  on  next  page. 
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Fig.  C.2(b)  Detailed  Comparison  of  Instantaneous  Vorticity  Contours  Computed  Using 
ARC2D  (Top)  and  ARC2D/8E  (Bottom),  tU/c  =  0.051,  Re  =  400,000.  Continued 
from  previous  page. 
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Comparison  of  Instantaneous  Vorticity  Contours  Computed  Using  ARC2D  (Top) 
and  ARC2D/8E  (Bottom),  tU/c  =  0.060,  Re  =  400/)00.  (Spacing  between  contours 
is  An  =  80;  dashed  lines  indicate  negative  levds). 
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Fig.  C.4  Time-Step  Sensitivity  of  Surface  Flow  Topology  Evolution  Calculated  with 
ARC2D,  and  Comparison  with  ARC2D/8E  Results.  (Lines  are  zero  skin  friction 
contours;  At^  is  time-step  based  on  freestream  speed  of  soimd). 


PS1  (ARC2D/8E) 


Fig.  C.5  Re5niolds  Number  Scaling  of  the  Time  of  Appearance  of  Incipient  Primary  (PSD 
and  Secondary  (SSI,  SS2,  SS3,  SS4)  Separation  Bubbles. 


Fig.  C6  Comparison  of  Vorticity  Fields  Near  the  Time  Beyond  the  End  of  Stage  II  for 
Re  =  50,000,  Re  =  100,000,  Re  =  200,000,  and  Re  =  400,000.  The  (Vortex-Induced) 
Formation  of  Eruptive  Plumes  of  Vorticity  is  Visible  at  the  three  highest 
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Fig.  C7  Instantaneous  Vortidty  Field  Beyond  Stage  n  at  tU/ c  =  0.078,  Re  =  400,000. 
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ieyond  Stage  n  at  tU/c  =  0.081,  Re  =  400,000. 
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FURTHER  RESULTS  ON  THE  REYNOLDS  NUMBER  SCALING 
OF  INCIPIENT  LEADING  EDGE  STALL 


In  Appendix  B  (AIAA  Paper  No.  94-2339),  four  Reynolds  numbers  were  investigated. 
These  are  (based  on  chord):  Re  =  50,000,  Re  =  100,000,  Re  =  200,000,  and  Re  =  400,000.  All 
simulations  were  for  a  fixed  angle  of  attack  a  =  15®,  and  a  freestream  Mach  number 
M„  =  0.2.  The  present  appendix  examines  further  results,  obtained  for  Re  =  800,000,  and 
how  these  affect  the  conclusions  drawn  in  Appendix  B.  These  results  will  be  presented  in 
AIAA  Paper  No.  95-0780. 

Grid  Convergence  Studies.  The  results  shown  here  are  computed  for  a  NACA0012  airfoil 
with  a  rounded  trailing  edge  (radius  of  curvature  r/c  =  0.0092).  This  ensures  that  the  CFL 
requirements  with  the  present  O-grid  are  not  dictated  by  the  trailing  edge  region,  which  is 
irrelevant  to  this  study.  The  scaling  relationships  discovered  in  Ref.  D1  provided  guidance 
in  the  design  of  the  computational  grids  at  Re  =800,000,  resulting  in  highly  clustered  meshes 
in  the  leading  edge  region  of  the  airfoil.  The  region  of  highest  circumferential  clustering 
was  located  between  x/c  =  0.002  and  x/c  =  0.025.  In  that  region,  the  circumferential  spacing 
was  imiform,  with  either  470, 235,  or  1 1 8  points,  depending  on  the  mesh. 

Three  computational  meshes  were  constructed  for  the  grid  convergence  studies 
conducted  at  Re  =  800,000.  Each  mesh  is  assembled  as  the  union  of  two  zones.  The  grid  in 
each  zone  is  designed  using  NASA  Ames's  GRIEXjEN  hyperbolic  grid  generation  program. 
This  particular  release  of  GRIDGEN  (Version  2.0)  allows  for  non-uniform  initial  spacing  in 
the  marching  direction  (i.e.,  normal  to  the  surface).  The  ability  to  prescribe  a  variable  initial 
grid  spacing  plays  an  essential  role  in  the  alleviation  of  CFL  requirements  at  the  trailing 
edge.  For  example,  in  all  of  the  grids  used  here,  the  initial  trailing  edge  spacing  normal  to 
the  surface,  Ayxg/  is  40  times  greater  than  the  initial  normal  grid  spacing  at  the  leading  edge 
(Ay^g/c  =  1.25x10'^).  An  example  of  a  computational  grid  used  in  this  study  is  given  in 
Fig.  D.l.  For  this  grid,  666  points  are  used  in  the  circumferential  direction,  470  of  which  are 
placed  within  2.3%  of  chord  (525  in  the  first  3%),  on  the  suction  surface.  A  detail  of  the  first 
(or  "inner")  grid  zone  is  illustrated  in  the  bottom  graph  of  Fig.  D.l.  The  inner  grid  zone 
(666x32  points)  surrounds  the  airfoil,  and  is  characterized  by  a  relatively  mild  stretching 
factor,  Ay|^.i/Ayj  =  1.065.  The  outer  grid  zone  extends  radially  approximately  seven  chord 
lengths  from  the  airfoil,  with  a  stretching  factor  of  1.17.  For  all  meshes  the  number  of  points 
in  the  direction  normal  to  the  airfoil  was  91,  with  an  initial  spacing  at  the  leading  edge 
AyLE  =  1.25x10"^.  The  value  of  Ay^g  was  conservatively  chosen  to  be  one  half  of  the  value 
used  at  Re  =  400,000  and  Re  =  200,000  (see  Ref.  Dl). 

It  is  worthwhile  noting  that,  for  the  present  problem,  normal  grid  spacing  is  not  an 
issue  and  can  be  based  on  standard  boundary  layer  scaling  arguments.  What  is  at  issue  is 
the  streamwise  resolution  needed,  because  of  the  possibility  of  self-focusing,  eruptive 
plumes  of  vorticity.*^^  The  effect  of  circumferential  grid  spacing  was  investigated  on  three 
computational  meshes  (666x91, 394x91,  and  258x91),  and  is  shown  in  Figs.  D.2,  D.3,  and  D.4. 
The  goal  of  this  exercise  is  to  find  out  when  local  grid  independence  is  attained.  By  "local" 
we  mean  "in  the  first  three  percent  of  chord,"  where  all  of  the  important  physics  of  incipient 
separation  take  place  (at  this  Reynolds  number  and  angle  of  attack). 

Skin  friction  was  used  as  a  sensitive  indicator  of  the  relative  accuracy  of  the  numerical 
solutions.  Figure  D.2  illustrates  the  convergence  of  skin  friction  distribution  between 
0.5%  chord  and  3.0%  chord,  as  the  circumferential  grid  spacing  is  reduced  from 
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Ax  in  =  2.5x10"^  to  AXjnin  =  6.25x10'^.  The  results  are  shown  at  four  instants  in  time  ranging 
from"  U/c  =  0.0375  to  tU/c  =  0.051.  The  two  finest  grids  yield  virtually  identical  results  up 
to  0.048  <  tU/c  <  0.051,  just  prior  to  the  onset  of  the  third  occurrence  of  a  secondary 
separation  bubble  (SS3).  As  previously  observed  at  Re  =  200,000,  the  accuracy  of  the 
numerical  solutions  is  subject  to  a  form  of  sensitivity  to  initial  conditions;  given  enough 
time,  two  solutions  will  always  diverge  from  one  another,  regardless  of  the  resolution  used. 
This  is  more  clearly  illustrated  in  Fig.  D.3,  where  the  evolution  of  the  skin  friction  error  is 
shown  as  a  function  of  time.  The  skin  friction  "error"  is  defined  as  the  L-1  norm  of  the 
difference  between  a  given  skin  friction  distribution  and  the  "converged"  skin  friction 
distribution  at  some  fixed  instant  in  time.  The  "converged"  Cf  distribution  is  taken  to 
correspond  to  the  highest  grid  resolution  case,  which  is  characterized  by 
Ax  in/c  =  6.25x10'^.  This  is,  of  course,  an  imperfect  assumption  but  nothing  substantially 
better  can  be  done.  It  is  important  to  realize  that  the  type  of  grid  convergence  study 
considered  here  is  local,  not  only  spatially  but  temporally.  The  spatial  and  temporal 
"regions"  where  the  solution  is  determined  to  be  accurate  are  the  factors  used  in  the  choice 
of  an  appropriate  computational  mesh.  Because  the  666x91  mesh  gave  grid-independent 
results  up  to  and  beyond  to  point  of  vortex  formation,  it  was  considered  adequate  for  the 
present  purpose,  i.e.;  to  obtain  grid-independent  results  for  incipient  leading  edge 
separation  at  Re  =  800,000. 

For  completeness,  the  scaling  of  the  skin  friction  error  as  a  function  of  circumferential 
grid  spacing  is  indicated  for  various  times  in  Fig.  D.4.  The  purpose  of  such  plots  is  to 
determine  the  effective  order  of  accuracy  and  degree  of  convergence  of  the  method.  The 
characteristic  slopes  associated  with  formal  first,  second,  third,  and  fourth-order  accuracy 
are  also  indicated  (dashed  lines)  for  reference.  At  small  times,  the  scaling  of  error  indicates 
an  effective  order  of  accuracy  which  lies  between  third-order  and  fourth-order.  Similar 
results  were  obtained  for  the  normal  grid  spacing  study  of  Ref.  D1  (Fig.  4,  Appendix  B)  at 
Re  =  200,000.  This  result  is  expected  because  initially  the  process  is  essentially  a  boundary 
layer  process  (i.e.,  viscous-dominated),  and  the  viscous  fluxes  in  the  code  are  computed  with 
fourth-order  accuracy.  In  general,  however,  lower  than  fourth-order  accuracy  is  expected 
because  the  formal  accuracy  of  the  differencing  scheme  degrades  as  one  approaches  the 
boundary.  The  effective  order  of  accuracy  on  the  skin  friction  remains  about  third-order  up 
to  tU/c  =  0.0435,  which  will  be  shown  to  be  beyond  the  time  of  primary  vortex  formation. 
By  tU/c  =  0.051,  the  accuracy  degrades  to  first  order.  This  is,  again,  an  indication  that  the 
results  should  be  interpreted  with  caution  for  tU/c  >  0.048. 

ArtificisI  Dissipation.  At  the  high  laminar  Reynolds  numbers  considered  here,  a  possible 
concern  is  the  extent  to  which  artificial  dissipation  alters  the  effective  viscosity  of  the  flow. 
Therefore,  the  ratio  between  the  added  flux  due  to  artificial  dissipation,  and  the  total 
viscous  flux,  Tu  (corresponding  to  the  "natural"  viscosity),  was  computed  at  every  point 
within  the  region  of  interest.  The  result  for  Re  =  800,000  is  shown  in  contour  form  in 
Fig.  D.5,  at  a  time  tU/c  =  0.06.  Typically,  the  ratio  Tad/'J'//  is  closely  tied  (a)  to  the  high- 
wavenumber  content  of  the  solution,  and  (b)  to  the  grid  spacing.  It  will  be  shown  that,  at 
tU /c  =  0.06,  multiple  fully-formed  vortices  are  present  in  the  solution  and  that  the  surface 
flow  topology  is  complex,  involving  at  least  four  distinct  regions  of  secondary  flow 
separation  and  two  regions  of  tertiary  separation.  Also,  the  reattachment  occurs  around 
x/c  =  0.025.  Therefore  the  "active"  vortical  region  of  the  flow  extends  roughly  between  1.0%, 
chord  and  2.5%  chord,  and  extends  less  than  0.5%  in  the  direction  normal  to  the  surface.  In 
this  region,  the  dissipation  flux  ratio,  0.01%.  The  results 

discussed  above  are  representative  of  the  overall  study.  Similar  tests  were  conducted  at  all 
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of  the  Reynolds  numbers  considered  here  (50,000  .  Re  .  800,000).  In  none  of  ‘^ese  tests  did 
^^er  exceed  0.1%  over  the  region  of  interest.  Thus,  there  is  a  good  degr^  o 
coiSidence  that  the  effective  Reynolds  number  is  indeed  800,000  for  the  results  presented  m 

this  section. 

Flowfield  nacrripfinn  at  Re  =  800.000.  A  time  sequence  of  the  developing  vorticity  field  is 
shown  in  Fig.  D.6.  At  the  very  early  stages  (tU/c  <  0.007),  a  Stokes  layer  orms  a  er 
imposing  the^no-slip  boundary  condition.  The  strong  adverse  pressure  gradient  eventually 
low  moLnlum  flU  near  the  wall  to  change  f  '“e  a  a 

nrlmarv  recircuUtion  region  (designated  "PSl")  is  initiated  around  tU/c  =  0-016,  at  a 
focation  indicated  in  Fig.  D.6(b).  This  is  foUowed  by  a  period  of  spatial  development  of  a 
thin  recirculation  bubble  near  the  surface,  the  beginning  and  end  of  which  are  dMignated  y 
ft” in  Figs.  D.6(c),  D.6(d),  and  D.6(e).  The  bubble  is  initially  rough  y 

symiJ^etric  in  shape  about  the  location  of  maximum  normal  displacement,  but  avenl^l  y 
bSlesTo pslded  as  vorticity  accumulates  in  the  reattachment  region.  This  is  visible,  for 

irrstance,  at  tU/c  =  0.0375  (end  of  Stage  I). 

The  evolution  of  the  flowfield  is  fairly  rapid  after  this  point.  Similarly  to  the  Reynolds 
number  cases  previously  analyzed,  a  primary  vortex  forms  in  the  downstream  region  of  the 
Sary  relation  zLe.  This  vorL,  in  turn,  induces  at  the  wall  a  region  of  seconda^ 
flow  separation.  The  extent  of  this  secondary  separation  is  indicated  by  the  s)^^ols  SSI  m 
Fies  D^6(f)  through  D.6(i).  As  time  progresses,  multiple  vortices  are  formed  as  a  result  o 
fXctfon  and  o^f  shear  layer  instability  and  roll-up.  Two  vortices  are  Present  ^ 
tU/c  =  0  045  three  vortices  at  tU/c  =  0.0525,  four  vortices  at  tU/c  =  0.06,  and  hve  vortices  at 
tU/c  -  0  0675  As  in  the  case  of  the  primary  vortex,  the  successive  shear  layer  vortices 
i!lc;  Lw  secondary  separation  regions.  These  multiple  secondary  recirculation  zones 
evenlually  sfreugften  to  the  point  that  they  too  become 

mav  induce  tertiary  separation.  Using  the  same  notation  as  in  Appendix  B,  the  success 
Tenons"  f  secondary^low  are  denoted  SSI,  SS2,  SS3,  etc.,  based  on  their^der  of 
appearance.  A  similar  convention  is  used  for  regions  of  tertiary  flow  s^^ration, 
aSd  so  on.  Note  that,  for  clarity,  contours  of  negative  vorticity  are  omitted  from  Fig.  D.6. 

Overall  the  flowfield  evolution  at  Re  =  800,000,  as  depicted  in  Figs.  D.6(a) 
through  D  6(h),  bears  strong  qualitative  similarity  with  previously  observed  results  at  lowe 
Revnolds  numbers  For  tU/c  <  0.06,  the  main  differences,  aside  from  the  usual  spatial  and 
temporal  scales  involved,  are  the  complexity  of  the  flowfield  and  in 

successive  shear  layer  vortices  present  in  the  solution  at  one  time.  For  tU/c  >  0.06,  severa 
qualitative  differences  were  observed,  as  compared  to  the  results  of  lower 
lumbers  At  Re  =  400,000,  some  differences  were  already  observed,  notably:  the  existence 

«fevLx-indu«dempHons(AppendixC).  However, these 

to  the  behavior  at  lower  Reynolds  numbers  occurred  after  the  rupture  of  the  FJ^ary 
feeding  sheet.  At  Re  =  800,000  further  differences  are  observed;  these  are  described  below. 

The  detailed  events  taking  place  between  tU/c  =  0.0675  (Fig.  D.6(i))  and  tU/c  =  0  0^ 
(Fig.  D.6(j))  are  depicted  in  Fig.  D.7,  using  a  smaller  time  increment 

on  the  experience  gained  at  the  lower  Reynolds  numbers,  the  appearance  of  the  vorticity 
contoursPn  Fig  D.6(i)  suggests  the  most  tenuous  link  in  the  vorticity  feeding  sheet  to  be 
between  the  sePond  and  third  vortices,  rather  than  upstream  of  the  primary  ' 

a  rapid  redistribution  of  the  vorticity  takes  place  between  the  two  vortices,  '^th  much  o 
PorLity  of  the  second  vortex  flowing  into  the  primary  vortex.  Simultaneously,  the 
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secondary  and  tertiary  flow  induced  by  the  primary  vortex  begin  to  lift  away  from  the 
surface,  "pinching  off  what  had  begun  to  look  like  the  merging  of  the  first  and  second 
primary  vortices.  This  is  clearly  visible  at  tU/c  =  0.069  and  tU/c  =  0.0698  in  Fig.  D.7.  This 
"pinching  off  halts  the  merging  process,  eventually  reversing  it  (i.e.,  splitting  the  primary 
structure)  at  tU/c  =  0.0705.  This  is,  of  course,  a  form  of  the  vortex-induced  eruptions 
previously  observed  at  Re  =  100,000,  Re  =  200,000,  and  Re  =  400,000.  There  is,  however,  a 
key  difference:  at  Re  =  800,000  this  eruption  occurs  underneath  a  vortex,  primarily,  rather 
than  between  two  vortices. 

For  tU/c  2  0.071,  the  primary  vortex  structure  is  no  longer  connected  to  its  feeding 
sheet  and  begins  to  convect  downstream.  Meanwhile,  the  split  upstream  portion  of  the 
merged  structure  (referred  to  as  "Vortex  2,"  for  convenience)  is  engulfed  into  the  third 
vortex.  Thus,  although  a  vortex-induced  eruption  is  also  observed  at  Re  =  800,000,  many  of 
the  details  of  the  interactions  surrounding  this  event  differ  from  the  lower  Reynolds  number 
cases.  In  particular,  the  key  event  is  the  vorticity  eruption  itself,  rather  than  the  rupture  of 
the  feeding  sheet.  For  50,000  <  Re  <  400,000  these  events  invariably  coincided.  In  contrast, 
at  Re  =  800,000,  a  rupture  of  the  vortex  sheet  between  vortices  does  not  occur.  Instead, 
vorticity  eruption  takes  place,  and  is  sufficiently  powerful  to  split  the  primary  structure  in 
two.  For  consistency  with  the  terminology  previously  used  in  this  work,  the  event  taking 
place  at  tU/ c  =  0.0705  (Fig.  D.7(c))  will  be  referred  to  as  the  point  of  "rupture,"  in  spite  of  the 
fact  that  its  significance  is  somewhat  different  than  at  lower  Reynolds  numbers.  The  results 
of  Fig.  D.7  support  previous  conjectures  that  the  vortex  sheet  rupture  observed  in  other 
calculations  is  indeed  triggered  by  the  eruptive  plume  of  vorticity. 

Dependence  of  Skin  Friction  on  Reynolds  Number.  Figures  D.8  through  D.12  depict  the 
space-time  evolution  of  the  skin  friction  at  Re  =  50,000,  Re  =  100,000,  Re  =  200,000, 
Re  =  400,000,  and  Re  =  800,000.  Figures  D.8  through  D.12  characterize  the  instantaneous 
surface  flow  topology,  which  typically  correlates  well  with  the  key  features  of  the  overall 
topological  structure  of  the  flow.  Results  on  the  pressure  are  presented  in  a  subsequent 
section  of  this  Appendix. 

The  qualitative  similarity  between  the  surface  flow  topologies  at  the  different  Reynolds 
numbers  can  be  seen.  As  in  Fig.  D.6,  the  flow  is  always  characterized  initially  by  the 
formation  and  development  of  a  primary  separation  bubble.  A  vortex  subsequently  forms 
in  the  downstream  region  of  the  bubble,  inducing  strong  negative  velocities  (large  negative 
skin  friction)  and  the  formation  of  a  secondary  separation  bubble.  This  scenario  is  very 
much  a  recursive  one:  in  time,  a  secondary  vortex  forms  in  the  "downstream"  region 
(upstream  with  respect  to  the  main  flow)  of  the  secondary  bubble,  inducing  strong  positive 
velocities  (large  positive  skin  friction)  and,  possibly,  the  formation  of  a  tertiary  separation 
bubble.  The  presence  of  successive  multiple  vortices  in  the  shear  layer  is  clearly  identified  at 
the  surface  by  their  induced  secondary  separation  regions,  i.e.,  the  number  of  regions  of 
positive  skin  friction  at  any  given  x/c.  As  previously  discussed,  the  number  of  formed 
shear  layer  vortices  is  seen  to  increase  with  Reynolds  number. 

Comparing  Figs.  D.8  through  D.12,  the  solution  that  was  carried  out  the  farthest  in 
time,  relative  to  flowfield  development,  is  at  Re  =  400,000.  There,  as  well  as  at  Re  =  200,000, 
a  break-away  point  exists  (the  rupturing  event)  after  which  the  primary  vortex  and 
subsequent  shear  layer  vortices  start  convecting  downstream  at  a  more  rapid  pace.  This  is 
indicated  by  the  change  in  slope  of  the  contours  on  the  x-t  diagrams.  In  Fig.  D.ll,  a  pairing 
can  be  identified  by  the  cross-over  of  negative  skin  friction  regions  (around  tU/c  =  0.05  and 
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x/c  =  0.095).  This  corresponds  to  the  event  depicted  in  Fig.  C.13  of  Appendix  C.  At 
Re  =  800,000  (Fig.  D.12),  the  large  negative  skin  friction  region  associated  with  the  primary 
vortex  appears  to  suddenly  end  at  tU/ c  =  0.073  and  x/c  =  0.025,  instead  of  convecting 
downstream,  as  in  the  lower  Reynolds  number  cases.  This  behavior  is  the  surface  signature 
of  the  event  depicted  in  Fig.  D.7(f),  namely  the  ejection  of  vorticity  and  the  lifting-away  of 
the  primary  vortex  structure.  Although  such  lift-away  is  also  noticed  at  Re  =  400,000 
(Fig.  D.ll),  it  does  not  occur  before  the  shedding  of  the  primary  vortex.  Thus,  the  behavior 
of  the  flowfield  at  Re  =  800,000  does  fundamentally  and  qualitatively  differ  from  the 
behavior  observed  for  Re  s  400,000.  However,  these  qualitative  differences  all  pertain  to  the 
later  stages  of  the  incipient  stall  process  namely,  after  the  formation  of  the  primary  stall 
vortex.  We  may  speculate  that,  as  the  Reynolds  number  is  increased  further,  the  complex 
vortex  system  will  stay  increasingly  confined  to  the  leading  edge  region  before  the  eruption 
occurs.  This  is  in  contrast  to  the  behavior  observed  for  Re  <  400,000,  in  which  the  vortex 
system  was  found  to  convect  first,  before  breaking  away  from  the  surface. 

The  space-time  skin  friction  diagrams  of  Figs.  D.8-D.12  can  be  simplified  to  some  extent 
by  considering  only  the  zero  contours,  corresponding  to  a  sign  reversal  of  the  skin  friction. 
In  this  manner,  the  lines  delimit  adjacent  regions  of  forward  and  reverse  flow  at  the  surface. 
This  permits  the  identification  of  primary,  secondary  and,  if  applicable,  tertiary  separation 
regions.  The  resulting  surface  flow  topology  diagrams  are  compared  in  Fig.  D.13  for  all  five 
Reynolds  numbers  investigated  in  this  study.  The  symbols  "PS,  SS,  and  TS  refer  to 
primary  separation,  secondary  separation,  and  tertiary  separation,  respectively.  The 
numeral  that  follows  each  of  these  symbols  indicates  the  temporal  order  of  appearance,  so 
that  for  example,  "SS2"  designates  the  second  occurrence  in  time  of  a  secondary  separation 
region,  etc.  The  surface  flow  topologies  at  the  various  Reynolds  numbers  are  seen  to  differ 
in  scale,  both  spatially  and  temporally.  There  is  also  a  difference  between  Re  =  800,000  and 
the  Re  s  400,000  cases:  at  lower  Reynolds  number  cases  are  characterized  by  a  bulk  motion 
of  the  leading  edge  vortex  system,  corresponding  to  convection;  in  contrast,  the  smface  flow 
patterns  at  Re  =  800,000  exhibit  less  bulk  downstream  motion,  and  more  of  a  "break,"  as  seen 
in  Fig.  D.12.  Apart  from  these  differences,  the  important  point  of  Fig.  D.13  is  the  remarkable 
degree  of  similarity,  including  many  of  the  details,  between  surface  flow  topologies  across 
the  Reynolds  number  range.  In  light  of  the  newest  results  at  Re  =  800,000,  the  qualitative 
simUarity  noted  in  Ref.  D1  is  found  to  hold,  except  after  the  point  of  vorticity  eruption.  This 
is,  again,  because  the  vortex-induced  eruption  at  high  Reynolds  number  does  no  longer 
coincide  with  the  rupture  of  the  vorticity  feeding  sheet  as  in  Appendix  B. 

Dependence  of  Pressure  Coefficient  on  Reynolds  Number.  The  space-time  evolution  of  the 
pressure  field  at  the  surface  is  documented  in  Figs.  D.14  through  D.18.  Close  examination  of 
the  pressure  coefficient  contours  reveals  many  of  the  same  conclusions  drawn  from  the  skin 
friction.  Low  pressure  areas  associated  with  the  presence  of  vortices  and  high-pressure 
ridges  associated  with  flow  impingement  immediately  downstream  of  these  vortices  are 
noticeable  features  of  the  flow  at  all  Reynolds  numbers.  The  relief  of  the  low  pressure 
region  downstream  of  the  leading  edge,  as  the  boundary  layer  and  separation  bubble  form, 
is  also  apparent. 

By  taking  successive  "horizontal  cuts"  of  Figs.  D.14-D.18,  one  obtains  a  more 
quantitative  description  of  the  pressure  distribution  and  its  temporal  evolution.  The. 
resulting  Cp  distributions  are  shown  in  Figs.  D.19-D.23.  To  further  illustrate  the  qualitative 
similarity  ot  the  flowfield  s,  the  total  time  period  shown  in  each  figure  corresponds  to  6.4% 
beyond  the  time  of  rupture.  For  each  Reynolds  number  this  period  is  divided  equally  into 
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17  time  increments.  In  the  first  frame  of  each  figure,  the  decrease  of  the  peak  suction  and 
the  relief  of  the  pressure  gradient  downstream  of  the  leading  edge  are  apparent.  The 
development  of  the  primary  separation  bubble  is  associated  with  a  lowering  of  the  pressure 
in  the  downstream  part,  followed  by  a  relative  increase  in  the  pressure  due  to  the  stagnation 
point.  The  associated  adverse  pressure  gradient  is  subsequently  characterized  by  a  sharp 
steepening  during  and  after  the  period  of  formation  of  the  primary  stall  vortex.  This  stage 
typically  corresponds  to  the  second  frame  in  each  figure.  Finally  (third  frame),  multiple 
primary  and  secondary  vortices  form,  greatly  complicating  the  interpretation  of  the  surface 
pressure  signature.  Eventually,  eruption  occurs  (corresponding  to  the  time-trace  before 
last),  and  the  vortex  system  starts  convecting  (last  trace  in  each  figure).  Again,  the 
qualitative  similarity  between  Reynolds  numbers  is  observed. 

We  conclude  this  section  with  the  calculation  of  the  streamwise  pressure  gradient.  For 
dynamic  stall,  the  pressure  gradient  is  by  far*^^  the  largest  contributor  to  the  production  of 
vorticity  at  the  wall.  In  the  absence  of  blowing /suction  or  surface  acceleration,  the  pressure 
gradient  is  directly  proportional  to  the  diffusive  normal  vorticity  flux  at  the  surface. It  is 
useful,  therefore,  to  attempt  to  obtain  insights  into  the  vorticity  production  process  from 
inferences  based  on  the  pressure  gradient.  This  approach  was  used,  e.g.,  by  Acharya  and 
Metwally  (Ref.  D5)  to  show  the  existence  of  a  "spike"  in  vorticity  flux  near  the  leading  edge 
of  a  pitching  airfoil  at  88,000  <  Re  <  120,000.  For  convenience,  the  streamwise  derivative  of 
the  Cp  distribution,  dCp/ax,  was  directly  computed,  rather  than  the  true  tangential 
derivative,  dCp/ds.  For  the  present  purpose,  this  makes  no  difference,  since  the  same  airfoil 
at  the  same  angle  of  attack  is  used  at  all  Reynolds  numbers.  The  result  is  shown  in 
Figs.  D.24-D.28.  The  complicated  spatio-temporal  patterns  which  follow  the  evolution  of 
the  vortices  are  the  most  noticeable  feature,  but  are  of  little  interest  in  terms  of  discerning 
the  origin  of  vortex  formation,  at  least  in  the  sense  argued  in  Ref.  D2. 

The  results  obtained  in  the  present  study  do  not  support  the  interpretation  that  the 
sudden  appearance  of  an  initial  eruptive  plume  of  vorticity  (which  would  presumably  be 
associated  with  a  spike  in  vorticity  flux)  is  responsible  for  the  formation  of  the  dynamic  stall 
vortex.  Instead,  the  present  problem  is  characterized  by  the  continuous  presence  of  a  strong 
vorticity  flux  near  the  leading  edge.  This  vorticity  flux  is  indeed  the  source  of  the  vorticity 
feeding  sheet.  With  time,  this  vorticity  sheet  lifts  from  the  surface,  creating  a  thin  separation 
bubble,  and  eventually  leads  to  the  formation  of  the  primary  stall  vortex.  In  this  process,  no 
qualitative  differences  were  found  across  the  investigated  Reynolds  number  range.  A 
consistent  behavior  for  all  Reynolds  numbers  is  the  initial  decrease  in  the  magnitude  of  the 
vorticity  flux  (Figs.  D24-D.281  This  initial  decrease  is  associated  with  the  formation  of  the 
viscous  boundary  layer,  which  changes  the  effective  shape  of  the  airfoil  and,  thus,  relieves 
the  pressure  gradient.  This  is  fundamentally  different  from  the  occurrence  in  time  of  a  spike 
in  the  vorticity  flux,  and  may  be  an  artifact  of  the  present  model  problem.  We  note, 
however,  that  for  Reynolds  numbers  exceeding  100,000  (Figs.  D.25-D.28),  the  leading  edge 
vorticity  flux  eventually  increases  again.  Whether  this  has  any  connection  with  the 
formation  of  a  spike  (in  the  sense  of  Ref.  D2)  is,  at  best,  speculative  at  this  time.  In  any 
event,  such  strengthening  takes  place  (for  the  present  problem)  after  the  formation  of  the 
dynamic  stall  vortex,  which  is  more  closely  connected  to  the  initial  flux  after  the  imposition 
of  the  no-slip  boundary  condition. 

Incipient  Vortex  Formation.  The  manner  in  which  the  point  of  vortex  formation  is 
determined  was  alluded  to  in  Appendix  B.  The  present  section  attempts  to  clarify  the 
methods  used  to  determine  both  incipient  vortex  formation  and  the  relation  between  the 
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topological  "center"  and  a  vortex,  based  on  the  results  at  Re  =  800,000.  "V ortex  formation" 
can  be  defined  as  the  transition  between  an  elongated  recirculation  bubble  (characterized 
primarily  by  counter-flowing  streams  separated  by  a  thin  shear  layer)  and  a  "fully-formed" 
vortex  structure  (characterized  by  significant  momentum  transfer  normal  to  the  wall  and 
the  existence  of  a  low-pressure  core). 

In  an  attempt  to  quantitatively  define  the  point  of  vortex  formation,  consider  the 
temporal  evolution  of  the  pressure.  Fig.  D.29.  At  short  times,  the  pressure  field  closely 
resembles  the  potential  flow  pressure  field,  with  the  exception  of  a  small  region  of  zero 
normal  pressure  gradient  near  the  wall  (the  boundary  layer).  This  can  be  seen,  for  example, 
in  Fig.  D.29(a).  At  tU/c  >  0.03  the  pressure  gradually  lowers  for  x/c  <  0.02  while 
simultaneously  increasing  downstream  of  that  point  (see  also  Fig.  D.23).  Around  a  time 
tU/c  =  0.037,  the  topology  of  the  pressure  field  undergoes  a  transformation,  namely  the 
existence  of  pressure  extrema  away  from  the  airfoil  surface.  The  closed  contours  of 
minimum  pressure  coincide  with  the  core  of  the  nascent  vortex.  At  that  point,  the  evolution 
of  the  flowfield  is  very  nonlinear:  the  low  pressure  core  contributes  to  the  further 
accumulation  of  vorticity  in  the  vortex.  The  increased  vorticity  in  turn  increases  the  radial 
pressure  gradient,  and  so  on,  until  the  vortex  is  fully  formed.  By  tU/c  =  0.045  (Fig.  D.29(g)), 
a  second  shear  layer  vortex  has  also  formed,  as  indicated  by  the  presence  of  low  pressure 
cores.  In  the  following,  the  transition  of  the  pressure  field  from  a  boundary-layer-like 
behavior  (i.e.,  dp/dy  =  0)  to  one  characterized  by  the  existence  of  closed  low-pressure 
contours  away  from  the  surface  is  defined  as  the  point  of  vortex  formation  (in  this  case: 
tvpU/c  =  0.037). 

An  important  element  is  to  establish  the  origin  of  the  vortex,  before  it  forms.  To  this 
end,  and  following  the  methodology  used  in  Refs.  D6  and  Dl,  we  use  Perry  and  Chong's^ 
critical  point  theory.  Knight  and  Choudhuri  (Ref.  D6)  have  shown  how  critical  point  theory 
can  be  used  as  a  sensitive  means  to  analyze  the  topology  of  the  velocity  field.  In  particular, 
these  authors  were  able  to  pinpoint  the  origin  of  the  dynamic  stall  vortex  for  a  pitching 
airfoil  at  Re  =  10,000.  Using  similar  concepts,  it  can  be  shown  that  the  critical  points  are 
"centers,"  characterized  by  essentially  zero  dilation,  D  =  du/ dx  +  dv/ dy,  and  positive 
Jacobian,  J  H  Ou/ax)(av/8y)  -  (au/dy)(dv/dx).  Fig.  D.30  exemplifies  the  method  used  to 
determine  the  location  of  the  topological  center.  In  this  example  (Re  =  400,000,  tU/c  =  0.03), 
there  is  no  vortex.  The  inset  in  Fig.  D.30  shows  the  actual  aspect  ratio  of  the  developed 
recirculation  bubble.  The  grey  contours  in  Fig.  D.30  (both  inset  and  main  figure)  are  low 
Mach  number  contours  (M  <  0.01)  which  help  demarcate  both  the  surface  flow  and  the  shear 
layer.  The  main  figure  is  obtained  by  rotating  and  stretching  the  coordinates,  as  indicated 
by  the  (x-y)  reference  arrows.  Also  present  in  the  figure  are  instantaneous  streamlines 
which  delineate  the  location  of  the  center. 


Close  examination  of  such  pictures  at  successive  instants  in  time  reveals  that  the 
presence  of  a  topological  center  emerges  very  shortly  after  the  onset  of  surface  flow  reversal. 
Immediately  after  applying  the  no-slip  condition,  the  Jacobian  is  zero  at  the  surface, 
indicating  a  "pure  shear"  situation.  After  the  onset  of  flow  reversal,  the  value  of  the  jacobian 
J  nears  zero  ("degenerate  center"),  then  rapidly  becomes  positive  ("center").  For  a  fairly  long 
period  of  time  corresponding  to  the  development  of  the  recirculation  bubble,  the  center 
moves  downstream,  within  the  primary  bubble,  but  does  not  correspond  to  a  vortex.  The 
later  stages  of  this  evolution  are  exhibited,  for  example,  in  Fig.  D.31,  which  shows  both 
instantaneous  pressure  contours  and  streamtraces  for  the  Re  =100,000  case.  As  previously 
noted,  the  presence  of  a  vortex  is  characterized  by  significant  momentum  transfer  normal  to 
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the  wall  and  the  existence  of  a  low  pressure  core.  In  Fig.  D.31,  this  occurs  between 
tU/c  =  0.054  and  tU/c  =  0.06.  In  all  of  the  cases  examined,  the  topological  center  is  clearly 
the  precursor  of  the  stall  vortex.  Its  trajectory  is  indicated  by  the  dark  symbols  in  the 
surface  topology  plot  of  Fig.  D.32.  The  behavior  at  Re  =  800,000  is,  again,  qualitatively 
similar  to  the  behavior  observed  at  lower  Reynolds  numbers,  for  instance:  Fig.  8  of 
Appendix  B,  for  Re  =  100,000.  In  particular,  until  the  point  of  vortex  formation  (marked 
"VF"  in  Fig.  D.32),  the  topological  center  moves  towards  the  two-thirds  downstream 
location  within  the  primary  bubble.  After  that  point,  its  trajectory  correlates  well  with  the 
downstream  edge  of  the  primary  and  induced  secondary  bubbles. 

The  coimection  between  flowfield  topology  and  vortex  formation  is  further  explored  in 
Fig.  D.33.  In  this  figure,  the  initial  temporal  increase  of  the  Jacobian  value  at  the  critical 
point  is  characteristic  of  the  emergence  of  a  center,  as  discussed  in  Appendix  B.  After  the 
point  of  vortex  formation  the  Jacobian  increases  further;  however,  this  growth  shows  no 
sign  of  correlation  with  the  onset  of  secondary  flow  features,  as  was  suggested  at 
Re  =  100,000  (Fig.  10,  Appendix  B).  The  dark  symbols  in  Fig.  D.33  show  how  vortex 
formation  is  characterized  by  a  drop  in  core  pressure.  This  feature  was  observed  at  all  of  the 
investigated  Reynolds  numbers. 

The  following  is  a  summary  of  the  topological  evolution  of  the  flow  at  Re  =  800,000. 
For  a  fairly  long  portion  of  time  during  Stage  I,  (tU/ c  <  0.038  in  Fig.  D.32),  the  vorticity  field 
topology  corresponds  to  a  lifted  vorticity  layer  beneath  which  a  thin  region  of  reverse  flow 
exists.  During  this  portion,  the  location  of  the  topological  center  "rides"  the  vorticity  line, 
gradually  moving  from  front  to  back  and  settling  at  about  two  thirds  of  the  bubble  extent. 
In  contrast  to  the  rather  benign  behavior  of  the  vorticity  during  that  portion,  the  end  of 
Stage  I  is  marked  by  a  sudden  change  from  boundary  layer-like  behavior  (negligible  normal 
pressure  gradient)  to  vortex-like  behavior,  characterized  by  the  presence  of  closed  pressure 
contours  associated  with  a  low  pressure  core  (Fig.  D.29).  At  later  stages  of  the  evolution  of 
the  flow,  local  pressure  and  vorticity  fields  act  to  reinforce  each  other.  The  vorticity 
accumulates  further  in  the  low  pressure  region,  and  the  increased  vorticity  results  in  a 
further  lowering  of  the  pressure.  This  interaction  results  in  a  formed  vortex  and  marks  the 
beginning  of  what  was  referred  to  earlier  as  Stage  II.  At  that  point,  vorticity  and  pressure 
fields  are  well-correlated. 

Reynolds  Number  Scalinp.  To  complete  the  Reynolds  number  scaling  investigation  initiated 
in  Ref.  Dl,  the  data  obtained  at  Re  =  800,000  was  processed  in  a  fashion  similar  to  that 
described  in  Appendix  B.  Figure  D.34  shows  the  scaling  with  respect  to  Reynolds  number 
of  the  time  instants  corresponding  to  the  onset  of  flow  reversal  (PSD,  vortex  formation,  and 
shear  layer  "rupture."  The  previously  obtained  power-law  curve  fits  for 
50,000  <  Re  <  400,000  are  indicated  by  the  lines  in  Fig.  D.34.  The  data  at  Re  =  800,000 
establish  a  noticeable  departure  from  previous  power-law  scalings.  In  particular,  both 
vortex  formation  and  rupture  occur  "late,"  relative  to  the  expected  scaling  based  on  the 
lower  Reynolds  number  data.  Of  course,  "rupture"  of  the  feeding  vorticity  sheet  per  se  did 
not  take  place  at  Re  =  800,000.  Instead,  the  vortex-induced  eruption  ended-up  splitting  the 
merged  primary  vortex.  This  event  is  depicted  in  Fig.  D.7  (frames  (b)  and  (c)).  It  is  not  clear 
at  this  point  whether  the  delay  of  the  rupturing  event  at  Re  =  800,000  can  be  accounted  for 
by  the  altered  definition  of  "rupture"  alone.  This  is,  in  fact,  unlikely,  because  of  the 
coincidence  of  the  feeding  sheet  rupture  with  vorticity  eruption  for  Re  <  400,000.  One  of  the 
key  differences  at  the  higher  Reynolds  number  is  the  presence  of  pairing/merging  events 
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before  shedding  occurs.  Additionally,  even  the  point  of  vortex  formation  appears  to  be 
somewhat  delayed  at  Re  =  800,000.  This  issue  is  explored  further  below. 

In  an  attempt  to  understand  further  the  scaling  properties  of  incipient  dynamic  stall 
with  respect  to  Reynolds  number,  the  boundary  layer  thickness  was  calculated.  Its  scaling  is 
shown  in  Fig.  D.35.  The  boundary  layer  thickness,  was  defined  as  follows.  The 
temporal  growth  of  the  boundary  layer  (initially  a  Stokes  layer)  was  monitored  at  the 
leading  edge,  and  was  recorded  once  it  reached  steady-state  (typically  for  tU/ c  <  0.006).  For 
convenience,  was  defined  as  the  normal  distance  from  the  wall  at  which  the  vorticity 
is  10%  of  its  value  at  the  wall.  Thus,  6^  le  is  not  the  vorticity  thickness  but,  rather,  a 
measure  of  boundary  layer  thickness  based  on  the  normal  damping  of  the  vorticity  profile. 
For  reference,  EE  corresponds  in  the  case  of  a  Blasius  boundary  layer  to 
ij  =  y(U„/vx)^/^  =  4.52,  or  6.8  momentum  thicknesses  ("U/U„"  =  0.98).  The  slope  of  the 
power-law  curve  fit  (solid  line  in  Fig.  D.35)  is  approximately  -1/2,  as  expected  for  a  laminar 
boimdary  layer. 

The  Reynolds  number  scaling  shown  in  Fig.  D.34  can,  thus,  be  re-examined  using 
different  length  scales  than  the  airfoil  chord.  Indeed,  the  chord  C  (including  in  the  definition 
of  the  Reynolds  number)  is  rather  irrelevant  to  the  processes  involved  in  incipient  leading 
edge  stall.  Three  alternative  length  scales  are  considered:  b  (the  leading  edge  radius  of 
curvature),  6^  ee  (upstream  boundary  layer  thickness),  and  the  viscous  penetration  length, 
(vt)l/2^  The  latter  is  a  natural  candidate  for  the  initial  stages  of  the  flow  development, 
which  are  essentially  similar  to  those  occurring  in  Stokes'  first  problem.  The  re-scaled  times 
tU/b,  tU/6f^EE'  tU/fvt)^'^^  are  shown  in  Fig.  D.36.  The  Reynolds  number  scaling  of  the 
various  times  of  onset  is,  of  course,  not  expected  to  change  by  using  b  as  the  length  scale, 
since  the  airfoil  is  the  same  and  since  b  does  not  change  with  Reynolds  number.  The  bottom 
two  frames  of  Fig.  D.36  suggest  the  following  observations:  (1)  temporal  normalization 
using  boundary  layer  (6^  le/U)  or  viscous  ((vt)'^^/U)  time  scales  does  not  collapse  the 
data,  and  (2)  the  relative  tardiness  of  the  vortex  formation  time  appears  to  be  eliminated  by 
using  either  (6^  ee/U)  or  ((vt)^/^/U).  Of  note  is  the  almost  perfect  fit  of  the  primary 
separation  time  of  onset,  when  non-dimensionalized  using  the  viscous  length  scale  (bottom 
frame).  The  results  of  Fig.  D.36  indicate  that  none  of  the  characteristic  length  scales 
considered  here  make  any  of  the  recorded  events  Reynolds  number  independent. 
Furthermore,  no  single  scale  can  be  expected  to  collapse  all  of  the  data.  This  points  to  the 
existence  of  a  hierarchy  of  scales,  a  view  somewhat  consistent  with  Ref.  D8. 

Finally,  using  the  additional  candidate  length  scales  (6^  ee  (vt)^^^),  we  re-examine 
the  scaling  of  the  times  of  onset  of  the  primary  and  secondary  recirculation  bubbles,  adding 
the  Re  =  800,000  case  to  the  data  presented  in  Appendix  B.  This  is  given  in  Fig.  D.37.  Again, 
the  data  appear  to  scale  according  to  simple  power-laws  of  the  Reynolds  number.  No  single 
time  scale  collapses  any  of  the  measured  variables.  With  the  exception  of  the  rupturing 
event  at  Re  =  800,000,  it  may  be  concluded  from  the  qualitative  self-similarity  exhibited  in 
Figs.  D.36  and  D.37  that  the  present  study  does  not  support  the  hypothesis  of  a  Reynolds 
number  bifurcation  under  laminar  flow  conditions. 

Auxiliary  Results.  For  completeness,  the  temporal  evolution  of  the  overall  force  coefficients 
for  lift,  drag,  and  moment  are  included  in  this  report  (Figs.  D.38-D.40).  Note  that  the  results 
in  Figs.  D.38-D.40  are  presented  for  tU/c  >  0.  No  attempt  is  made  to  interpret  these  results 
in  connection  with  the  leading  edge  dynamic  stall  phenomena  described  in  Appendices  B, 
C,  and  D.  This  is  because  Ce,  Ce),  and  C,„  are,  to  a  large  extent,  determined  by  conditions  at 
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the  trailing  edge.  Since  the  present  geometry  considers  a  rounded  trailing  edge,  the 
unsteady  Kutta  condition  allows  for  the  rear  stagnation  point  to  move,  both  in  time  and 
with  Reynolds  number.  The  details  of  this  region  of  the  flow,  although  critical  for  lift  and 
drag,  are  largely  irrelevant  in  the  context  of  the  present  study.  A  natural  question  is  to  ask 
whether  the  choice  of  an  airfoil  geometry  was  the  correct  one,  considering  the  specific  goals 
of  the  research.  With  the  presently  clustered  computational  grids  (notably  at  Re  =  400,000 
and  800,000),  the  cost  of  the  computational  grid  points  associated  with  the  airfoil  (O-grid) 
geometry  accounts  only  for  a  fraction  of  the  total  computational  cost.  This  fraction 
decreases  with  increasing  Reynolds  number.  In  addition,  the  use  of  an  airfoil  ensures  that 
the  proper  pressure  gradient  is  imposed,  without  having  to  experiment  with  approximate 
boundary  conditions. 

Finally,  we  conclude  this  appendix  with  the  documentation  of  the  temporal  evolution 
of  the  vorticity  fluxes  at  the  leading  edge.  Although  for  different  conditions  (Re  =  12,000 
pitching  airfoil).  Appendix  A  addressed  the  issue  of  predictability  of  the  leading  edge 
vorticity  fluxes  before  and  during  the  dynamic  stall  process,  by  means  of  a  modified  version 
of  indicial  theory.  One  of  the  important  results  was  the  conclusion  that,  although  difficulties 
remained  for  predicting  the  fluxes  directly,  it  appeared  that  a  prediction  of  the  vorticity 
accumulation  based  on  indicial  theoretical  concepts  is  possible.  For  nonlinear  mohons 
(Aa  -  5*),  this  prediction  is  a  relative  one,  permitting  accurate  prediction  of  one  airfoil 
maneuver  from  the  knowledge  of  another.  The  potential  applications  of  these  findings  were 
never  fully  exploited  in  the  course  of  the  present  study,  in  part  because  of  changes  in  A.R.O. 
research  priorities  part-way  through  the  research.  In  an  attempt  to  begin  bridging  this  gap, 
the  last  three  figures  document  the  evolution  of  the  integrated  tangential  and  normal  fluxes 
of  vorticity  in  the  present  semi-impulsive  problem  at  Re  =  800,0(X). 

Figure  D.41  depicts  the  temporal  evolution  of  the  normally-integrated  tangential 
vorticity  flux,  4>t,  at  0.0032  <  x/c  <  0.0282.  At  all  locations,  the  flux  is  strongly  dominat^  by 
the  initial  conditions  of  the  present  problem;  notably,  all  fluxes  initially  decrease  until  the 
point  of  vortex  formation.  Figure  D.42  documents  the  temporal  evolution  of  the  normal 
flux  of  vorticity,  «!►„,  tangentially  integrated  between  x/c  =  0.0  and  x/c  =  0.03.  The  normal 
flux  is  shown  at  two  locations:  at  the  wall  (solid  line)  and  at  0.03%  chord  away  from  the 
wall.  As  in  the  case  of  <I>t,  the  behavior  of  the  normal  flux  at  the  wall  is  initially  dominated 
by  its  relaxation  from  the  initial  conditions.  The  temporal  growth  of  the  boundary  layer  can 
also  be  correlated  to  the  initial  increase  (0  <  tU/c  <  0.006)  in  <!>„  away  from  the  wall  This 
flux  peaks  and  subsequently  decreases  just  prior  to  the  point  of  vortex  formation  ("VF'). 
Figure  D.43  shows  the  evolution  of  the  vorticity  accumulation  in  the  first  3%  of  chord.  Also 
indicated  for  reference  are  the  times  of  onset  corresponding  to  primary  separation  (PSl), 
vortex  formation  (VF),  secondary  separation  (SSI,  SS2,  SS3,  SS4),  and  rupture  (R).  After  its 
initial  decrease,  it  is  noteworthy  that  the  vorticity  accumulation  begins  to  increase  around 
tU/c  =  0.03,  which  happens  to  coincide  with  the  simultaneous  strengthening  of  the 
topological  center  and  the  reduction  in  core  pressure,  as  shown  in  Fig.  D.33.  More  complete 
analyses  of  pitching/oscillating  airfoil  data  need  to  be  performed  before  one  can  establish  a 
correlation  between  vorticity  accumulation  and  incipient  vortex  formation.  Should  such  a 
correlation  exist,  then  the  results  of  Appendix  A  imply  that  a  semi-analytical  prediction  of 
dynamic  stall  is  possible. 
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error 


Fig.  D.3  Temporal  Evolution  of  Skin  Friction  Error  With  Respect  to  Finest  Grid  Solution 
(Re  =  800/)00). 
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Fig.  D.4  Scaling  of  Skin  Friction  Error  as  a  Function  of  Circumferential  Grid  Spacing  at 
Four  Instants  in  Time  (Re  =  800,000).  (Dashed  lines  indicate  slopes  associated 
with  formal  first,  second,  third,  and  fourth  order  accuracy). 


ratio 


Temporal  Evolution  of  the  Vorticity  Field  at  Re  =  800.000.  (Arrows  indicate 
locations  of  skin  friction  sign  reversal).  (Continued  on  next  page). 
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Fig.  D.8  Contour  Plot  Illustrating  the  Space-Time  Evolution  of  Skin  at  Re  =  50,000.  (Black 
lines  indicate  zero  skin  friction). 


Re  =  100,000 


(Black  lines  indicate  zero  skin  friction). 
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(Black  lines  indicate  zero  skin  fiiction). 


Re  =  400,000 


Contour  Plot  Illustrating  the  Space-Time  Evolution  of  Skin 
(Black  lines  indicate  zero  skin  friction). 


Re  =  800,000 
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(Black  lines  indicate  zero  skin  friction). 
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Fig.  D.13  Zero  Skin  Friction  Contours  Illustrating  the  Similarity  of  Surface  Flow 
Topologies  Across  the  Reynolds  Number  Range. 


Re  =  50,000 


Fig.  D.14  Contour  Plot  Illustrating  the  Space-Time  Evolution  of  the  Pressure  Coefficient  at 
the  Surface^  Re  =  50,000. 


Re  =  100,000 


the  Surface,  Re  =  lOO/XX). 


Re  =  200,000 


Fig.  D.16  Contour  Plot  Illustrating  the  Space-Time  Evolution  of  the  Pressure  Coefficient  at 
the  Surface,  Re  =  200/)00. 


Re  =  400,000 


Fig.  D.17  Contour  Plot  Illustrating  the  Space-Time  Evolution  of  the  Pressure  Coefficient  at 
the  Surface,  Re  =  400,000. 


Re  =  800,000 
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Fig.  D.18  Contour  Plot  Illustrating  the  Space-Time  Evolution  of  the  Pressure  Coefficient 
the  Surface,  Re  =  800,000. 
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Fig.  D.19  Pressure  Coefficient  Distribution  as  a  Function  of  Time,  Re  =  50,000. 
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Fig.  D.20  Pressiire  Coefficient  EHstribution  as  a  Function  of  Hme,  Re  =  100,000, 
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Fig.  D.22  Pressure  Coeffident  Efetribution  as  a  Function  of  Time,  Re  =  400,000. 


Re  =  800,000 


Fig.  D.23  Pressure  Coeffident  Distribution  as  a  Function  of  Time,  Re  -  800,000. 
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Fig.  D.26  Contour  Plot  Illustrating  the  Space-Time  Evolution  of  the  Streamwise  Pressure 
Gradient  at  the  Surface,  Re  =  200,000. 
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Gradient  at  the  Surface,  Re  =  400,000. 
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Fig.  D.29  Temporal  Evolution  of  the  Pressure  Field 
Illustrating  Vortex  Formation  at  Re  =  800,000. 
(6(p/pJ  indicates  the  value  spacing  between 
successive  contour  lines). 
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30  Schematic  Illustration  of  the  Methodology  Used  to  Determine  the  Location  of  the 
Topological  Center  at  tU/c  =  0.03,  Re  =  400,000.  (The  inset  depicts  the  extent  of 
the  un-rotated,  un-scaled  recirculation  region  using  iso-contours  the  Mach 
number,  for  M  s  0.01.  Lines  in  the  main  graph  are  instantaneous  streamtraces). 


Temporal  Evolution  of  the  Pressure  Field  and  of  Streamtraces  Illustrating  the 
Connection  between  the  Topological  Center  and  the  Formation  of  the  Primary 


SIDGIiH 


Vortex  Formation 


Temporal  Evolution  of  Key  Features  Used  to  Describe  the  Primary  Vortex 
History  at  Re  =  800,000.  (O  Jacobian,  ■  pressure  at  the  core). 


Scaling  of  the  Upstream  Boundary  Layer  Thickness  with  Reynolds  Number. 
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Fig.  D.36  Alternate  Reynolds  Number  Scaling  of  Onset  Times  for  Primary  Separation, 
Vortex  Formation,  and  Rupture,  Using  Various  Characteristic  Length  ^ales 
(Top;  Radius  of  Curvature,  Center  Boundary  Layer  Thickness,  Bottom:  Viscous 
Length  Scale). 


Fig.  D.37  Alternate  Reynolds  Number  Scaling  of  Onset  Times  for  Incipient  Primary  (PSD 
and  Secondary  (SSI,  SS2,  SS3,  SS4)  Separation  Bubbles,  Using  Various 
Characteristic  Length  Scales  ((a)  Airfoil  Chord,  (b)  Radius  of  Curvature, 
(c)  Boimdary  Layer  Thickness,  (d)  Viscous  Length  Scale). 
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Fig.  D.38  Lift  Coefficient  Time  Histories. 
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Fig.  D.39  Drag  Coefficient  Time  Histories. 


Fig.  D.40  Pitching  Moment  Coefficient  Time  Histories. 
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♦  V 


Time  History  of  the  Integrated  Normal  Vorticity  Flux  at  the  Leading  Edge/ 
Re  =  800,000. 


